Code archives/Algorithms/Handy Exotic Maths

This code has been declared by its author to be Public Domain code.

Download source code

Handy Exotic Maths by Streaksy2010
I think a lot of this hasn't been done before but is still quite useful.

I made all of these functions but the following exceptions:
Direction#()
LinesCross()
InTriangle()

If anyone knows the authors of those functions (I may have rewrote parts of them), please comment. So sorry I couldn't do it myself. I hate it when people aren't creditted for the contributions.

Anyway... Have a look through and see if any of it appeals to you. Cheers =D
Global sm$=Chr$(34) ;convenience speech-mark character (quote)

Global floatyrez=10000
Global TriCenterX2D#,TriCenterY2D#

Global RotatedPointX#,RotatedPointY#




; MATHS

;returns value between V1 and V2 with T as tween (0-1)
Function Between#(v1#,v2#,t#):dif#=v2-v1
Return v1+(dif*t)
End Function

;distance between 2 2D points
Function Distance#(x1#,y1#,x2#,y2#)
Return Sqr((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2))
End Function

;distance between 2 3D points
Function Distance3D#(x1#,y1#,z1#,x2#,y2#,z2#)
Return Sqr#((x1-x2)^2+(y1-y2)^2+(z1-z2)^2)
End Function

;returns the tween (0-1) of where between L and R the value M is
Function WhereBetween#(l#,r#,m#)
r=r-l:m=m-l
Return m/r
End Function

;returns the absolute difference between two values
Function Difference#(val1#,val2#)
Return Abs(val1-val2)
End Function

;returns whether an integer is an even number
Function Even(v)
If Abs(v)=1 Then Return 0
Return Not (Abs(v) And 1)
End Function

;returns whether an integer is an odd number
Function Odd(x)
Return Abs(x) And 1
End Function

;returns v1 to the power of v2.  (blitz's `^' operation is buggy with high values, and sometimes in general.  this was just as fast in the stress-test anyway)
Function power(v1,v2)
ov1=v1
If v2=0 Then Return 1
If v2=1 Then Return v1
For t=1 To v2-1
v1=v1*ov1
Next
Return v1
End Function

;rounds a float to the closest integer
Function Round(val#) 
dec#=val-Floor(val):If dec<.5 Return Floor(val) Else Return Ceil(val)
End Function


;curve a value based on its range and a curve type and amplitude (see my Curve() upload in the code archives for better example)
Function Curve#(val#,min#,max#,typ=3,amp#=1)
val=val-min
max=max-min
If amp<>1 Then olval#=val
tween#=((val/max)*90)
If typ=<1 Then cos1#=Cos(tween-90):val=cos1*max					;smooth out
If typ=2 Then cos1#=1-Cos((tween)):val=cos1*max					;smooth in
If typ=3 Then cos1#=Cos(tween-90)*Sin(tween):val=cos1*max		;smooth in and out
If amp<>1 Then dif#=olval-val:val#=olval-(dif*amp) ;amplify
Return val+min
End Function

;returns the difference between two angles in degrees
Function AngleDifference#(angle1#,angle2#)
Return ((angle2 - angle1) Mod 360 + 540) Mod 360 - 180
End Function

;returns an angle after turning it to another angle a given amount (not as straight-forward as it sounds)
Function TurnTo#(ang#,angto#,turnamt#)
a1#=ang
a2#=ang
For t=1 To 180
a1#=a1#-1
a2#=a2#+1
If a1<0 Then a1=a1+360
If a2=>360 Then a2=a2-360
If Int(a1)=Int(angto) Then
ang=ang-turnamt
EndIf
If Int(a2)=Int(angto) Then
ang=ang+turnamt
EndIf
Next
Return ang
End Function

;returns the X vector(?) of an angle
Function AngleX#(aa#)
Return Cos(aa-90)
End Function

;returns the Y vector(?) of an angle
Function AngleY#(aa#)
Return Sin(aa-90)
End Function

;returns the direction that point B is from point A (I think I borrowed this from someone else... sincere apologies!!! but info should be spread)
Function Direction#(EnX#,EnY#,OtX#,OtY#)
If OtX#>EnX# And EnY#>=OtY# Then Return ATan((OtX#-EnX#)/(EnY#-OtY#))
If OtX#>=EnX# And OtY#>EnY# Then Return 90+ATan((OtY#-EnY#)/(OtX#-EnX#))
If EnX#>OtX# And OtY#>=EnY# Then Return 180+ATan((EnX#-OtX#)/(OtY#-EnY#))
If EnX#>=OtX# And EnY#>OtY# Then Return 270+ATan((EnY#-OtY#)/(EnX#-OtX#))
End Function

;limit a value to a minimum and maximum
Function Clamp#(V#,l#,u#)
Return Lim(v,l,u)
End Function

;limit a value to a minimum and maximum
Function Lim#(vl#,lw#,up#) ; A.K.A - Clamp
If vl<lw Then Return lw
If vl>up Then Return up
Return vl:End Function

;return a lower resolution float (reduce it to 2 decimal places)
Function Dec2#(v#)
Return Float(Int(v*10))/10
End Function

;return a lower resolution float (reduce it to 3 decimal places)
Function Dec3#(v#)
Return Float(Int(v*100))/100
End Function

;return a lower resolution float (reduce it to 4 decimal places)
Function Dec4#(v#)
Return Float(Int(v*1000))/1000
End Function

;returns a string of an integer with the given number of digits (places any zero on the front that it needs to)
Function PadNum$(val,digits=2)
v$=val:If Len(v$)<digits Then Repeat:v$="0"+v$:Until Len(v$)=digits
Return v$:End Function

;wrap a value between a minimum and a maximum (i love this.. never seen it done before :P )
Function Wrap(val,min,max)
If val=>min And val=<max Then Return val
If val>max Then val=val-min:max=max-min:max=max+1:Return min+(val Mod max)
If val<min Then min2=-min:max2=-max:max=min2:min=max2:val=-val:val=val-min:max=max-min:max=max+1:Return 0-(min+(val Mod max))
End Function

;wrap a float value... slower and clunkier but nessecities must...
Function WrapFloat#(val#,min#,max#)
span#=Abs(max-min)
.gain
cosh=cosh+1:If cosh>100000 Then RuntimeError "Wrap trap!"
If val>max Then val=val-span:Goto gain
If val<min Then val=val+span:Goto gain
Return val
End Function




;turn an integer into a string, using however many chracters (bytes) it takes
Function Word$(val)
Repeat
res$=res$+Chr(val Mod 255)
val=Int(val/255)
Until val=0
Return res$
End Function

;undo the affects of Word$(), turning a string back into an integer
Function UnWord(w$)
power=1
For t=1 To Len(w)
val=val+(Asc(Mid(w,t,1))*power)
power=power*255
Next
Return val
End Function





;turn a float into a string, using however many chracters (bytes) it takes.  The resolution (decimal places) is set at the top of the code
Function Floaty$(v#)
vv=Floor(v)
v2#=v-vv
v2s$=Str(v2):v2s=Right(v2s,Len(v2s)-2):If Len(v2s)>floatyrez Then v2s=Left(v2s,floatyrez)
vs$=Str(v)
w1$=xword(vv)
If v2<>0 Then w2$=xword(Int(v2s))
Return Chr(Len(w1))+w1+w2
End Function

;undo the affects of Floaty$(), turning a string back into a float
Function Unfloaty#(v$)
w=Asc(Left(Str(v),1))
v=Right(v,Len(v)-1)
w1$=xunword(Left(v,w))
If w<Len(v) Then w1=w1+"."+xunword(Right(v,Len(v)-w))
vv#=Float(w1)
Return vv#
End Function



Function xWord$(val) ;for floaty()
Repeat
res$=res$+Chr(val Mod 254)
val=Int(val/254)
Until val=0
res=Replace(res," ",Chr(255))
Return res$
End Function

Function xUnWord(w$) ;for unfloaty()
w=Replace(w,Chr(255)," ")
power=1
For t=1 To Len(w)
val=val+(Asc(Mid(w,t,1))*power)
power=power*254
Next
Return val
End Function


Function RectsOverlap%(x1, y1, w1, h1, x2, y2, w2, h2)
If x1>(x2+w2) Or (x1+w1)<x2 Then Return False
If y1>(y2+h2) Or (y1+h1)<y2 Then Return False
Return True
End Function





;I love this too.  also never seen it done.
;It returns a random value, but a float!  You can specify the resolution... so if resolution is .1 then it will return to 1 decimal place
Function Randm#(min=0,max,res#=1) ;res is resolution.  so random(0,1,.1) will be 0,.1,.2,.3,.4,.5,.6,.7,.8 or .9... random(5,20,5) will be 5,10,15,20
nums=(max-min)/res
Return (Rand(0,nums)*res)+min
End Function


;returns the angle after it has bounced off a "wall".  entrot# is the angle the "ball" is travelling, and Wallface# is the angle the wall is facing
Function Bounce2D#(entrot#,wallface#)
Return (wallface*2)-entrot
End Function


;get the center (really just a point inside) of a triangle (a bit rough but seems solid enough.  could use some improvement)
;it doesnt really return the centre (my trig is terrible.. someone help :P ) but it returns a point inside the triangle
;Resulting values are recorded in the global variables: TriCenterX2D# and TriCenterY2D#
Function GetTriangleCenter2D(x1#,y1#,   x2#,y2#,   x3#,y3#)
cx#=(x2+x3)/2
cy#=(y2+y3)/2
dir#=direction(x1,y1,cx,cy)
dis#=distance(x1,y1,cx,cy)
TriCenterX2D#=x1+  (      (   anglex(dir) * (dis/2)   )     )
TriCenterY2D#=y1+  (      (   angley(dir) * (dis/2)   )     )
End Function


;another one im rather proud of.  it returns the distance from a point to an infinite line.
Function DistanceFromLine#(x1#,y1#,x2#,y2#,   x#,y#)	;THE LINE IS EXTENDED FOREVER!  THIS JUST RETURNS A PARALLEL Y DIFFERENCE !
dir#=direction(x1,y1,x2,y2)+90
cx#=between(x1,x2,.5):cy#=between(y1,y2,.5)
pointdir#=direction(cx,cy,x,y)
pointdis#=distance(cx,cy,x,y)
newang#=pointdir-dir+180
Return (angley(newang)*pointdis)
End Function


;another favourite.  this rotates the new position of point A after it has been rotated around point B with the given angle.
;Resulting values are recorded in the global variables: RotatedPointX# and RotatedPointY#
Function RotatePoint(x#,y#,  ang#,  pivx#=0,pivy#=0)
dis#=distance(pivx,pivy,x,y)
dir#=direction(pivx,pivy,x,y)+ang
RotatedPointX#=pivx+(anglex(dir)*dis)
RotatedPointY#=pivy+(angley(dir)*dis)
End Function


;returns whether two lines cross.  this is one of the functions i found in the code archives.
;apologies to author for not noting who you are when i took this
Function LinesCross(x0#,y0#, x1#,y1#, x2#,y2#, x3#,y3# ,parralelcheck=1,pointtouch=0)   ;parralel means the lines dont cross but "overlap"
If pointtouch Then
If Abs(x0-x2)<.001 And Abs(y0-y2)<.001 Then Return True
If Abs(x1-x3)<.001 And Abs(y1-y3)<.001 Then Return True
EndIf
n# = (y0#-y2#)*(x3#-x2#) - (x0#-x2#)*(y3#-y2#)
d# = (x1#-x0#)*(y3#-y2#) - (y1#-y0#)*(x3#-x2#)
If Abs(d#) < 0.0001 
; Lines are parallel!
	If parralelcheck Then
	dir#=-direction(x0,y0,x1,y1)-90
	cx#=between(x0,x1,.5)
	cy#=between(y0,y1,.5)
	RotatePoint(x0,y0,dir,cx,cy)
	xx1#=RotatedPointX
	yy1#=RotatedPointY
	RotatePoint(x1,y1,dir,cx,cy)
	xx2#=RotatedPointX
	yy2#=RotatedPointY
	RotatePoint(x2,y2,dir,cx,cy)
	xx3#=RotatedPointX
	yy3#=RotatedPointY
	RotatePoint(x3,y3,dir,cx,cy)
	xx4#=RotatedPointX
	yy4#=RotatedPointY
	If xx1>xx2 Then xxxx=xx1:xx1=xx2:xx2=xxxx
	If xx3>xx4 Then xxxx=xx3:xx3=xx4:xx4=xxxx
	xx2=xx2-10
	xx3=xx3+10
	If xx2>xx3 And xx1<xx3 Then Return True
	If xx1<xx4 And xx2>xx4 Then Return True
Else
Return False
EndIf
Else
; Lines might cross!
Sn# = (y0#-y2#)*(x1#-x0#) - (x0#-x2#)*(y1#-y0#)
AB# = n# / d#
If AB#>0.0 And AB#<1.0
CD# = Sn# / d#
If CD#>0.0 And CD#<1.0
; Intersection Point
X# = x0# + AB#*(x1#-x0#)
Y# = y0# + AB#*(y1#-y0#)
Return True
End If
End If
; Lines didn't cross, because the intersection was beyond the end points of the lines
EndIf
; Lines do not cross!
Return False
End Function


;returns whether a point is within a given triangle (looks like i borrowed this too but i needed it for more functions below :D )
;wish i recorded the author of these functions so I could credit them...  so sorry... if anyone knows...
Function InTriangle(x0#,y0#   ,x1#,y1#,x2#,y2#,x3#,y3#)
b0# =  (x2 - x1) * (y3 - y1) - (x3 - x1) * (y2 - y1)
b1# = ((x2 - x0) * (y3 - y0) - (x3 - x0) * (y2 - y0)) / b0 
b2# = ((x3 - x0) * (y1 - y0) - (x1 - x0) * (y3 - y0)) / b0
b3# = ((x1 - x0) * (y2 - y0) - (x2 - x0) * (y1 - y0)) / b0 
If b1>0 And b2>0 And b3>0 Then Return True Else Return False
End Function

;return whether a line enters a triangle.  its probably not 100% perfect but it worked out great for a certain project
Function LineInTriangle(lx1#,ly1#,lx2#,ly2#,   tx1#,ty1#,tx2#,ty2#,tx3#,ty3#)
If intriangle(lx1,ly1,  tx1,ty1,tx2,ty2,tx3,ty3) Then Return True
If intriangle(lx2,ly2,  tx1,ty1,tx2,ty2,tx3,ty3) Then Return True
If linescross(lx1,ly1,lx2,ly2,   tx1,ty1,tx2,ty2,0,0) Then Return True
If linescross(lx1,ly1,lx2,ly2,   tx2,ty2,tx3,ty3,0,0) Then Return True
If linescross(lx1,ly1,lx2,ly2,   tx1,ty1,tx3,ty3,0,0) Then Return True
dis#=distance(lx1,ly1,lx2,ly2)
steps=5
For t#=0 To steps-1 Step 1
tw#=t/Float(steps-1)
x#=between(lx1,lx2,tw)
y#=between(ly1,ly2,tw)
If Abs(distancefromline(tx1,ty1,tx2,ty2,   x,y))>.00001 Then	;make sure its not too close to a triangle edge!
If Abs(distancefromline(tx2,ty2,tx3,ty3,   x,y))>.00001 Then
If Abs(distancefromline(tx1,ty1,tx3,ty3,   x,y))>.00001 Then
If intriangle(x,y,  tx1#,ty1#,tx2#,ty2#,tx3#,ty3#) Then Return True
EndIf
EndIf
EndIf
Next
End Function

;same as above but a bit less sophisticated.  I cant remember why I kept this or even what the difference is
Function LineInTriangleSimple(lx1#,ly1#,lx2#,ly2#,   tx1#,ty1#,tx2#,ty2#,tx3#,ty3#)
If intriangle(lx1,ly1,  tx1,ty1,tx2,ty2,tx3,ty3) Then Return True
If intriangle(lx2,ly2,  tx1,ty1,tx2,ty2,tx3,ty3) Then Return True
If linescross(lx1,ly1,lx2,ly2,  tx1,ty1,tx2,ty2,0,0) Then Return True
If linescross(lx1,ly1,lx2,ly2,  tx2,ty2,tx3,ty3,0,0) Then Return True
If linescross(lx1,ly1,lx2,ly2,  tx1,ty1,tx3,ty3,0,0) Then Return True
dis#=distance(lx1,ly1,lx2,ly2)
For t#=dis*.1 To dis*.9 Step .1
tw#=t/dis
x#=between(lx1,lx2,tw)
y#=between(ly1,ly2,tw)
If intriangle(x,y,tx1#,ty1#,tx2#,ty2#,tx3#,ty3#) Then Return True
Next
End Function

;i guess this is a more accurate version of LineInTriangle.  I made these a while ago so I forget :P
Function LineInTriangleHiDef(lx1#,ly1#,lx2#,ly2#,   tx1#,ty1#,tx2#,ty2#,tx3#,ty3#)
dis#=distance(lx1,ly1,lx2,ly2)
;dir#=direction(lx1,ly1,lx2,ly2)
If intriangle(lx2,ly2,  tx1#,ty1#,tx2#,ty2#,tx3#,ty3#) Then Return True
For t#=0 To dis Step .1
tw#=t/dis
x#=between(lx1,lx2,tw)
y#=between(ly1,ly2,tw)
If Abs(distancefromline(tx1,ty1,tx2,ty2,   x,y))>.00001 Then	;make sure its not too close to a triangle edge!
If Abs(distancefromline(tx2,ty2,tx3,ty3,   x,y))>.00001 Then
If Abs(distancefromline(tx1,ty1,tx3,ty3,   x,y))>.00001 Then
If intriangle(x,y,  tx1#,ty1#,tx2#,ty2#,tx3#,ty3#) Then Return True
EndIf
EndIf
EndIf
Next
End Function

Comments

Streaksy2010
I think I borrowed AngleDifference#() by the look of it. Doesn't look like something I could have made.


Warpy2010
Good stuff! Would benefit from some indentation, though.

My comments (there are lots, but you posted lots!):

WhereBetween needs a check for when r=l, and it should return r in that case.

Your Even and Odd functions use logical And instead of bitwise '&'. You should either change to (v & 1) or (v mod 2). The Abs() is unnecessary.

In the Curve function, note that Cos(a-90) = Sin(a), so you can write that instead. It's up to you though.

Here's a simpler TurnTo function:
Function TurnTo#(ang#,angto#,amount#)
	d# = AngleDifference(angto,ang)
	If Abs(d)<amount	'distance to 'angto' is less than 'amount'
		Return angto
	Else
		Return ang+Sgn(d)*amount
	EndIf
End Function


AngleX and AngleY don't need the '-90' bit, if you think of angle 0 as pointing to the right. Never mind this if you like to have 0 pointing upwards. The argument for zero pointing right is that it's a bit easier to think about as screen co-ordinates point right and downwards.

The Atan2 function can do exactly what your Direction function does.

Another way of doing 'Clamp' is:
Return Min(Max(vl,lw),up)


Your Dec2/Dec3/Dec4 functions all need an extra zero each.

Here's a simpler WrapFloat
Function WrapFloat#(val#,vmin#,vmax#)
	d  = 1+vmax - vmin
	a = (val-vmin)/d
	If val<vmin a=a-1
	Return val-a*d
End Function


For the triangle centre function, be aware that there are loads of different centres. The 'centroid' is the most intuitive one, and there are formulas for that all over the net.

To get a random float, use the Rnd function instead of Rand.

Here's a version of RotatePoint that uses a rotation matrix - basically the same as yours, but with fewer calculations
Function RotatePoint(x#,y#,  ang#,  pivx#=0,pivy#=0)
	x=x-pivx
	y=y-pivy
	RotatedPointX# = pivx+Cos(ang)*x-Sin(ang)*y
	RotatedPointY# = pivy+Sin(ang)*x+Cos(ang)*y
End Function





And that's all I've got to say about that. Don't be discouraged by my many nitpickings, I'm always in favour of more maths!


Wayne2010
what about :
Shortest distance between two lines in 3d ?
Line intersects cube in 3d?


Streaksy2010
Warpy, thanks so much =D I'm not the best mathmatician in the world, but I know what's useful. If there are better ways to do something then I want 'em. I'll look up that centroid thing soon too.

And Wayne, I have no clue how to do that stuff :P


Streaksy2010
Oh and Warpy... The custom random float thing is because I want to specify the resolution.

Hmm looks like I forgot to implement that :P


Streaksy2010
Here's a centroid function I just made. It was kinda obvious when I saw the diagram. :P Thanks again, Warpy.

THIS NEEDS TO BE ADDED TO LINESCROSS() AT THE PART IT FINDS THE INTERSECTION POINT:
LinesIntersectX=x
LinesIntersectY=y

LinesIntersectX and LinesIntersectY will have to be made global, too.


Function GetTriangleCenter2D(x1#,y1#, x2#,y2#, x3#,y3#)
midx2#=between(x2,x3,.5)
midy2#=between(y2,y3,.5)
midx3#=between(x3,x1,.5)
midy3#=between(y3,y1,.5)
cr=linescross(x1,y1,midx2,midy2 , x2,y2,midx3,midy3)
TriCenterX2D#=LinesIntersectX
TriCenterY2D#=LinesIntersectY
End Function


Floyd2010
The easy way to compute the centroid of a triangle is to take the "average" of the three corners, i.e.

Centroid_X = ( x1 + x2 + x3 ) / 3

and likewise for Y ( and Z if in 3D ).


Streaksy2010
Thats what I was doing I think, but that isn't a true centroid. If one corner is far away from the others, the other two will weigh the center too far to their side.


Floyd2010
There is only one centroid ( center of mass ). For 2D this is the point on which the object would balance.

For a triangle we pick any corner and draw the line segment to the midpoint of the opposite side. This segment is called a median. The triangle would balance if placed so this median was on a taut wire. That is intuitively clear if you think of the triangle as a collection of line segments parallel to the side. The median bisects each of them. The centroid must thus be on the median. The same argument applies to any vertex. So the centroid is the intersection point of the medians.

Finding this intersection point is the brute force way to get the centroid. But the "average vertex" method gives the same result.

As long as I'm rambling about centroids, here's another bit of trivia. On each median the centroid is exactly one third of the way from the side to the opposite vertex.


Streaksy2010
I did that median thing you mentioned 3 comments ago :P scroll up :P the corner-to-median intersection point is my new method. thats what you meant right?

but you get the same result from averaging the 3 coordinate sets? strange...


Floyd2010
I did that median thing you mentioned 3 comments ago

I know.

The average of three thing is often given as an exercise when you first about vectors.

First, take any fourth point and draw the three vectors from there to the three vertices. The average of these three vectors points to a location inside the triangle.

Next pick any vertex, construct the median ( using vector operations ) and find the position on it which is a third of the way from the opposite side to the vertex. Finally show that the average from the first step points exactly at this location.

Since this works for any choice of vertex what looked like three locations is really only one, and it must be the centroid since it is on all three medians.

As an alternative to all this theory you can simply try your method and the average of three method on a few test cases. Check that they give the same result.


Streaksy2010
Did a comparison...

Time taken to do line-intersect method 1000000 times: 171ms
Time taken to do half-median method 1000000 times: 271ms
Just taking averages of coordinates 1000000 times: 43ms

And, omg, yep the result is the same when you use the averages method. My bad :P





Function GetTriangleCenter22D(x1#,y1#, x2#,y2#, x3#,y3#) ;find intersection of 2 median lines
midx2#=(x2+x3)/2:midy2#=(y2+y3)/2
midx3#=(x3+x1)/2:midy3#=(y3+y1)/2
cr=linescross(x1,y1,midx2,midy2 , x2,y2,midx3,midy3)
TriCenterX2D#=LinesIntersectX
TriCenterY2D#=LinesIntersectY
End Function


Function GetTriangleCenter32D(x1#,y1#, x2#,y2#, x3#,y3#) ;find halfway-point on a median line
cx#=(x2+x3)/2
cy#=(y2+y3)/2
dir#=direction(x1,y1,cx,cy)
dis#=distance(x1,y1,cx,cy)
TriCenterX2D#=x1+ ( ( anglex(dir) * (dis/2) ) )
TriCenterY2D#=y1+ ( ( angley(dir) * (dis/2) ) )
End Function

Function GetTriangleCenter2D(x1#,y1#, x2#,y2#, x3#,y3#) ;average the 3 coordinates
TriCenterX2D#=(x1+x2+x3)/3
TriCenterY2D#=(y1+y2+y3)/3
End Function


_PJ_2010
The Power() function didn;t handle negative indices.
Function power#(v1,v2)
	ov1=v1
	If (Not(v2)) Then Return 1
	If v2=1 Then Return v1
	For t=1 To Abs(v2)-1
		v1=v1*ov1
	Next
	If (v2<0) Then v1=(1.0/v1)
	Return v1
End Function


Understandably, Fractional indeces are a whole different kettle of fish.


Code Archives Forum