For the graphics gurus, perhaps <smile>

General FreeBASIC programming questions.
Post Reply
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

This may be of no interest to the forum graphics gurus but since I don't know diddly-squat about graphics then there is no way I can tell. <smile>

I have put this in the General forum because it is not platform specific; a departure for me as most of my stuff is Windows based.

Imagine we have an object in two dimensional space moving in a straight line and we wish it to make a random turn. Clearly we can turn left from 0 to -180° degrees or turn tight from 0 to 180°. It may be that we wish to restrict a turn to ± 90° as otherwise it would appear unnatural.

Now we can use RND on the range [-90, 90] but any value in that range will have an equal chance of being chosen because RND is a uniform distribution. A spontaneous change in direction by, say, 90° may also be regarded as appearing unnatural.

Another approach is to use a different probability density function. We could use the normal distribution where a change in direction by, say, 90° would be highly unlikely, 80° more likely, 70° more likely still and so on on until we get to the greatest likely change in the range [-10°, 10°].

I am sure that many of you are familiar with the normal distribution and here it is. We have different shapes according to differing means and variances. The peak of the bell curve will cause a slight change in direction and the tails will see a greater change in direction but a very much less likely one.

Image

Better still, to my mind, is a variation of the normal distribution called the von Mises distribution. instead of two variables dictating the shape of the bell curve we have just one, referred to as kappa. Here it is.

Image

The tails of the normal distribution are asymptotic with the x axis and get increasing closer to the x axis as we move further away from x = 0. With von Mises the tails are bounded in the range [-π, π] or [-180°, 180°]. This suits our purpose.

What further interests me is what we get by differing values of kappa. With kappa equal to zero we get a uniform distribution, which is not what we want but as kappa increases we get the formation of a bell curve. When we get to kappa = 4 we are still outside the range [-90°, 90°] for some values, but not many. With kappa = 8 it looks we are now inside the range [-90°, 90°]. Beyond kappa = 8 our range tightens even more.

In the von Mises function below, because it is an approximation, we cannot use kappa = 0 as we get a division by zero and, as we all know, computers do not like that. In the following table I look at a variety of values of kappa starting with kappa = 0.0003.

10,000 calls to the function were made and counts taken in 10° steps from -90° to 90°. A count was also taken with < -90° and > 90°. Here is a table with some kappas.

Code: Select all

               kappa
 
   0.0003    4      8     16    24
 
    2492     38      0      0     0
-90
     282     47      0      0     0
-80
     273     53      5      0     0
-70
     287    131     27      3     0
-60
     288    256     68      2     0
-50
     263    411    205     21    10
-40
     289    668    453    189    66
-30
     271    905    926    628   402
-20
     266   1141   1510   1562  1519
-10
     258   1353   1839   2543  3072
  0
     283   1283   1815   2633  3025
 10
     305   1118   1452   1543  1486
 20
     272    977    882    689   372
 30
     278    679    509    153    46
 40
     290    424    204     28     2
 50
     252    223     70      6     0
 60
     272    125     26      0     0
 70
     270     81      5      0     0
 80
     272     42      3      0     0
 90
    2537     45      1      0     0
With a kappa = 0 we would get about 2500 < -90° and about 2500 > 90°. With kappa = 0.0003 I got 2492 and 2537 respectively. With kappa = 0 we would get about 278 in the 10° sub-sections. With kappa = 0.0003 we can see the uniform distribution pattern forming.

With kappa = 8 we only have one value from the 10,000 calls outside the range [-90°, 90°].

With kappa = 24 our range has decreased to [-50°, 50°]. Our change in direction of our object is now restricted to that range with the upper ends being the less likely to be returned.

The von Mises function uses some heavy arithmetic and the time spent in the function varies and is subject to three values of RND being used.

I did 20 tests of one million calls to establish the slowest frequency of returning values and got between 18.78KHz to 244.14KHz. These values are small compared with the frequency of a single RND but 18.78KHz is a return in 53 micro-seconds so may not have a bearing on a graphics program.

The function returns values in degrees but, obviously, that can be changed to your requirements.

The code is an adaptation of code by Michael Bell at the PowerBASIC forums.

Code: Select all

#define Pi 3.141592654
#define TwoPi 6.283185307
 
Declare Function vonMises( As Single ) As Single
 
Function vonMises( kappa As Single ) As Single
 
  '-----------------------------------------------------------------------------
  ' Returns random deviate in degrees from von Mises distribution with mean zero and a given kappa
  ' Uses Best & Fisher (1979) Method (Applied Statistics 28 pp.152-157)
  '-----------------------------------------------------------------------------
  Static as single tau, rho, r
  Static BeenHere As Long
  dim as double z, u1, u2, u3, f, c, theta
 
  If BeenHere = 0 Then
    ' initial quantities
    tau = 1+Sqr( 1+4*kappa*kappa )
    rho = ( tau-Sqr( 2*tau ) )/( 2*kappa )
    r = ( 1+rho*rho )/( 2*rho )
    BeenHere = -1
  End If
 
step1:
  u1 = RND : u2 = RND : u3 = RND
  z = COS( Pi*u1 )
  f = ( 1+r*z )/( r+z )
  c = kappa*( r-f )
 
step2:
  IF c*( 2-c ) - u2 > 0 THEN GOTO step4
 
step3:
  IF LOG( c/u2 ) + 1 - c < 0 THEN GOTO step1
 
step4:
  theta = SGN( u3-0.5 )*acos( f )
 
  IF theta > Pi THEN
    theta = theta - TwoPi
  ELSEIF theta <- Pi THEN
    theta = theta + TwoPI
  END IF
 
  Return theta*180/Pi
 
End Function
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

So, how accuarte is the vonMises function.

Code: Select all

For j = 1 to 10
  For i = 1 to 1000000
    total += vonMises(0.0003)
  Next
  Print total/1000000
  total = 0
Next
Theoretically, the average should zero.

The worst case I got was 0.133023109375 and since that is parts of 360° that translates to 0.00037 in parts of 1 so the arithmetic is quite good. Not as accurate as RND but in the same ball park.
paul doe
Moderator
Posts: 1733
Joined: Jul 25, 2017 17:22
Location: Argentina

Re: For the graphics gurus, perhaps <smile>

Post by paul doe »

deltarho[1859] wrote:This may be of no interest to the forum graphics gurus but since I don't know diddly-squat about graphics then there is no way I can tell. <smile>

I have put this in the General forum because it is not platform specific; a departure for me as most of my stuff is Windows based.
Well put, this is useful for several applications, not only graphics. I can see several uses for it, thanks. Mind if I steal it? ;)
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

paul doe wrote:Mind if I steal it? ;)
Of course I don't, provided you let me know what other uses you may have for it without giving away any trade secrets. <smile>

From a graphics perspective each moving object could have it's own kappa value. Of course, the graphics fan's imagination will be much greater than mine.
paul doe
Moderator
Posts: 1733
Joined: Jul 25, 2017 17:22
Location: Argentina

Re: For the graphics gurus, perhaps <smile>

Post by paul doe »

deltarho[1859] wrote:
paul doe wrote:Mind if I steal it? ;)
Of course I don't, provided you let me know what other uses you may have for it without giving away any trade secrets. <smile>
Hahaha you're a sharp one, indeed 8D

Love to discuss this later, as I'm busy right now ;)
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

If we change the vonMises return value to '(1 + theta/Pi)/2' we get a return range of [0,1] close to RND's [0, 1). For kappas > 0 there will be a bias to 0.5 and that bias will increase as kappa increases.

We could still be in the domain of graphics but with, effectively, the same range as RND then our imagination may get triggered.
dodicat
Posts: 7983
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: For the graphics gurus, perhaps <smile>

Post by dodicat »

This is Von's missus's postulate.
She didn't trust Von.
for e^-(x^2) she prefered LolRemez.
viewtopic.php?f=7&t=25897
A very shaky postulate though.

Code: Select all

Function f(x As Double) As Double
    ' approximating exp(-x*x)  [-5 To 5]
    Return _ 
    ((((((((((((((((((((((+1.474437987935631e-020)*x+1.939212847374133e-011)*x-1.993834680721019e-018)*x _ 
    -2.722158269845897e-009)*x _ 
    +1.158365607581166e-016)*x _ 
    +1.660219956710741e-007)*x _ 
    -3.774823092703865e-015)*x _ 
    -5.772195583993584e-006)*x _ 
    +7.551910689439813e-014)*x _ 
    +0.0001262198627201913)*x _ 
    -9.530345912825716e-013)*x _ 
    -0.001807738572678891)*x _ 
    +7.488807595474571e-012)*x _ 
    +0.01713941924820732)*x _ 
    -3.466122961144952e-011)*x _ 
    -0.1063740133180027)*x _ 
    +8.328600590366501e-011)*x _ 
    +0.4171114267672486)*x _ 
    -7.688585743517694e-011)*x _ 
    -0.9523132162264163)*x _ 
    +9.125150492293233e-012)*x _ 
    +0.9951565491299325) _ 
    
End Function

Dim As Long limit=200000
Redim As Single xp(1 To limit),yp(1 To limit),clr(1 To limit)
#define range(f,l) (Rnd*((l)-(f))+(f))

Screen 20
Window (-5,1)-(5,0)

For n As Long=1 To limit
    Var x=range(-5,5)
    xp(n)=Rnd*x
    yp(n)=f(x)
    clr(n)=1+Rnd*16
Next

Dim As Single L,M,t
Do
    Screenlock
    If Len(Inkey) Then m-=.1:t=Timer
    l+=.1
    Cls
    For n As Long=1 To limit
        yp(n)+=m*yp(n)
        Var xx= xp(n)+1*5*Sin(L)*yp(n)*Sin(L/30)
        Pset(xx,yp(n)),clr(n)
    Next
    Screenunlock
    Sleep 1,1
Loop Until m Andalso Timer-t>2

 
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

The things you graphics guys get up to - I love the collapse on clicking the close button.

That is what I meant by "the graphics fan's imagination will be much greater than mine."
DamageX
Posts: 130
Joined: Nov 21, 2009 8:42

Re: For the graphics gurus, perhaps <smile>

Post by DamageX »

Another approach is to use a different probability density function. We could use the normal distribution where a change in direction by, say, 90° would be highly unlikely, 80° more likely, 70° more likely still and so on on until we get to the greatest likely change in the range [-10°, 10°].
What if you get a random number between -1,1 and cube it. The range would remain -1,1 but the values would tend to be closer to 0. Then scale to +/-90 degrees.

Code: Select all

dim as double irnd,turn
irnd=RND*2-1
turn=irnd*irnd*irnd*90
deltarho[1859]
Posts: 4310
Joined: Jan 02, 2017 0:34
Location: UK
Contact:

Re: For the graphics gurus, perhaps <smile>

Post by deltarho[1859] »

Good idea. The thing about the normal distribution is that it is a natural phenomenon. Our brains may interpret a cubic approach as unnatural but I am only guessing.
Post Reply