|
atan2
Sept 23, 2022 6:55:03 GMT -5
Post by tenochtitlanuk on Sept 23, 2022 6:55:03 GMT -5
Many of us have found atan2, which takes variables y and x and gives an angle, to be very useful. I know I've coded it myself, and seen/used versions by Chris Iverson /Anatoly /B+ /Steelweaver /Rod /etc. Trouble is it is only needed infrequently. Recently, playing with arcs, I realised there are two common algorithm /definition versions in circulation. Checking on wikipedia- atan2 I find there is little consistency across programming languages. Up to 8 versions can be considered to exist, depending on whether you regard angles as being measured clockwise or anti-clockwise, and from which axis/direction regarded as zero. And you can define atan2 as giving results in the range 0 to 2 pi, or from -pi/2 to +pi/2!! No wonder I've sometimes got weird results and had to do a bit of fiddling! Just another symptom of my aging brain!
|
|
|
atan2
Sept 23, 2022 8:15:46 GMT -5
Post by Rod on Sept 23, 2022 8:15:46 GMT -5
Hmm... Here is a statement I would like to understand better.
In my ball collision code I use atan2 but according to the above statement I dont need to! So how to remove atan2 from this 2D collision code.
'------------------------------------------ '2D BILLIARD STYLE COLLISION '------------------------------------------ 'converted from a Blitz BASIC demo 'By Joseph 'Phish' Humfrey ' 'converted with the help of 'Anatoly and Stefan
nomainwin WindowWidth = 800+themeWidth WindowHeight = 600+themeHeight UpperLeftX = (DisplayWidth-WindowWidth)/2 UpperLeftY = (DisplayHeight-WindowHeight)/2 graphicbox #1.g, 0,0,800,600 open "2D Bumping" for graphics_nf_nsb as #1 #1 "trapclose [quit]"
#1.g "down ; fill black ; flush"
[setup] width=800 height=600 friction=.99 maxballs=15 dim ball(maxballs,6) 'x y delta x delta y radius mass b.x=1 'ball x b.y=2 'ball y b.dx=3 'ball delta x b.dy=4 'ball delta y b.r=5 'ball radius b.m=6 'ball mass
b.i=0 'ball index
'set out red balls in a triangle ballTriangleSize=5 for xloop = ballTriangleSize to 1 step -1 for yloop = 1 to xloop b.i=b.i+1 ball(b.i,b.x) = (5-xloop)*27 + 200 ball(b.i,b.y) = yloop*31-(xloop*31)/2.0 + 300 ball(b.i,b.dx)= 0 ball(b.i,b.dy)= 0 ball(b.i,b.r) = 15 ball(b.i,b.m) = 50 next next
'set out white cue ball cue=0 ball(cue,b.x) = 600 ball(cue,b.y) = 300 +20 ball(cue,b.dx)= -20 ball(cue,b.dy)= 2-int(rnd(0)*4+.5) ball(cue,b.r) = 15 ball(cue,b.m) = 50
#1.g "when characterInput [setup]" #1.g "setfocus"
timer 17,[updateballs] wait
[updateballs] for b.i= 0 to maxballs
'update positions ball(b.i,b.x)=ball(b.i,b.x)+ball(b.i,b.dx) ball(b.i,b.y)=ball(b.i,b.y)+ball(b.i,b.dy)
'gradually slow down ball(b.i,b.dx)=ball(b.i,b.dx)*friction ball(b.i,b.dy)=ball(b.i,b.dy)*friction
'------------------------------------------ 'COLLISION CHECKING '------------------------------------------ ' Check each ball against every other
for b= b.i to maxballs
' not against itself if b<>b.i then
collisionDistance = ball(b.i,b.r)+ball(b,b.r) actualDistance = sqr((ball(b,b.x)-ball(b.i,b.x))^2 + (ball(b,b.y)-ball(b.i,b.y))^2)
' the collision distance is the some of the radius of both balls ' if the actual distance appart is less than that then they have collided If actualDistance<collisionDistance Then
' calculate the angle the colliding ball is travelling at collNormalAngle=atan2(ball(b,b.y)-ball(b.i,b.y),ball(b,b.x)-ball(b.i,b.x))
' re position the balls as exactly touching, no intersection ' this prevents overlapping balls being seen as still colliding moveDist1=(collisionDistance-actualDistance)*(ball(b,b.m)/(ball(b.i,b.m)+ball(b,b.m))) moveDist2=(collisionDistance-actualDistance)*(ball(b.i,b.m)/(ball(b.i,b.m)+ball(b,b.m))) ball(b.i,b.x)=ball(b.i,b.x) + moveDist1*Cos(collNormalAngle+3.14159) ball(b.i,b.y)=ball(b.i,b.y) + moveDist1*Sin(collNormalAngle+3.14159) ball(b,b.x)=ball(b,b.x) + moveDist2*Cos(collNormalAngle) ball(b,b.y)=ball(b,b.y) + moveDist2*Sin(collNormalAngle)
'------------------------------------------ 'COLLISION RESPONSE '------------------------------------------ 'n = vector connecting the centers of the circles. 'we are finding the components of the normalised vector n nX=Cos(collNormalAngle) nY=Sin(collNormalAngle)
'now find the length of the components of each movement vectors 'along n, by using dot product. a1 = ball(b.i,b.dx)*nX + ball(b.i,b.dy)*nY a2 = ball(b,b.dx)*nX + ball(b,b.dy)*nY
'optimisedP = 2(a1 - a2) ' ---------- ' m1 + m2 optimisedP = (2.0 * (a1-a2)) / (ball(b.i,b.m) + ball(b,b.m))
'now find out the resultant vectors 'r1 = c1\v - optimisedP * mass2 * n ball(b.i,b.dx) = ball(b.i,b.dx) - (optimisedP*ball(b,b.m)*nX) ball(b.i,b.dy) = ball(b.i,b.dy) - (optimisedP*ball(b,b.m)*nY)
'r2 = c2\v - optimisedP * mass1 * n ball(b,b.dx) = ball(b,b.dx) + (optimisedP*ball(b.i,b.m)*nX) ball(b,b.dy) = ball(b,b.dy) + (optimisedP*ball(b.i,b.m)*nY)
end If end if next
'Simple Bouncing off walls. if ball(b.i,b.x)<ball(b.i,b.r) then ball(b.i,b.x)=ball(b.i,b.r) ball(b.i,b.dx)=ball(b.i,b.dx)*-0.9 end if if ball(b.i,b.x)>width-ball(b.i,b.r) then ball(b.i,b.x)=width-ball(b.i,b.r) ball(b.i,b.dx)=ball(b.i,b.dx)*-0.9 end if if ball(b.i,b.y)<ball(b.i,b.r) then ball(b.i,b.y)=ball(b.i,b.r) ball(b.i,b.dy)=ball(b.i,b.dy)*-0.9 end if if ball(b.i,b.y)>height-ball(b.i,b.r) then ball(b.i,b.y)=height-ball(b.i,b.r) ball(b.i,b.dy)=ball(b.i,b.dy)*-0.9 end if
next
[drawballs] #1.g "discard ; redraw" for b.i= 0 to maxballs if b.i=0 then #1.g "color white ; backcolor white" else #1.g "color red ; backcolor red" end if #1.g "place ";ball(b.i,b.x);" ";ball(b.i,b.y) #1.g "circlefilled "; ball(b.i,b.r) #1.g "color white ;\";b.i next
wait
[quit] close #1 end
function atan2(y,x) pi = asn(1) * 2 if x <> 0 then arctan = atn(y/x)
select case case x > 0 atan2 = arctan
case y>=0 and x<0 atan2 = pi + arctan
case y<0 and x<0 atan2 = arctan - pi
case y>0 and x=0 atan2 = pi / 2
case y<0 and x=0 atan2 = pi / -2 end select end function
|
|
|
atan2
Sept 23, 2022 8:34:16 GMT -5
Post by tsh73 on Sept 23, 2022 8:34:16 GMT -5
I fail to see how
A = Math.acos ( dot (v1, v2)/ (v1.length ()*v2.length ()) ); is any better then
collNormalAngle=atan2(ball(b,b.y)-ball(b.i,b.y),ball(b,b.x)-ball(b.i,b.x)) We don't have "dot()" or "length()" functions in LB so we will end up using more user defined functions Now, what for?
|
|
|
atan2
Sept 23, 2022 8:51:52 GMT -5
via mobile
Post by Rod on Sept 23, 2022 8:51:52 GMT -5
Well my trig is poor and I think I understand vectors and vector maths a little better. Maybe I should put the time in to better understand atan2.
But the statement did say the vector maths was quicker.
|
|
bplus
Full Member
Posts: 127
|
atan2
Sept 24, 2022 11:57:44 GMT -5
Post by bplus on Sept 24, 2022 11:57:44 GMT -5
Ha! My trig is far better than my vectors.
I base my x, y coordinates around point cx, cy like this: x = cx + r * cos(RadianAngle) y = cy + r * sin(RadianAngle) r = radial distance between x, y and cx, cy or r = SQR((x-cx)^2 + (y-cy)^2)
When RadianAngle = 0 then Cos(0) = 1, Sin(0) = 0 this graphs due East of point cx, cy When RadianAngle = Pi/2 = 90 degrees Cos(RadianAngle) = 0, Sin(RadianAngle) = 1 that graphs due South of cx, cy and moving from East to South is Clockwise around cx, cy and the angle increases from 0 to pi/2 (or if you prefer, 0 to 90 Degrees) an onwards around the point in 2*pi Radians or 360 Degrees.
For Atan2, I've used Andy Amaya's code with great results:
Function Atan2(y, x) 'Atan2 is a function which determines the angle between points 'x1, y1 and x2, y2. The angle returned is in radians 'The angle returned is always in the range of '-PI to PI radians (-180 to 180 degrees) '============================================================== 'NOTE the position of Y and X arguments 'This keeps Atan2 function same as other language versions '============================================================== If x = 0 Then If y < 0 Then Atan2 = -1.5707963267948967 Else Atan2 = 1.5707963267948967 End If Else chk = atn(y/x) If x < 0 Then If y < 0 Then chk = chk - 3.1415926535897932 Else chk = chk + 3.1415926535897932 End If End If Atan2 = chk End If 'thanks Andy Amaya End Function
Sometimes when comparing angles, you add 2*Pi to Radian Angle if < 0 so both angles compared are positive. Whether you use all pos or positive and negative angles, they graph to same point using the Trig functions I went over above.
And John is right, arc drawing in QB64 is backwards, so I had to create my own arc drawing Sub to draw arcs that match Basic's Trig functions. And I use this code in JB when translating:
global PI, DEG, RAD ' trig constants PI = acs(-1) DEG = 180 / PI RAD = PI / 180
'draws an arc with center at xCenter, yCenter, radius from center is arcRadius sub arc xCenter, yCenter, arcRadius, dAStart, dAMeasure 'notes: 'you may want to adjust size and color for line drawing 'using angle measures in degrees to match Just Basic ways with pie and piefilled 'this sub assumes drawing in a CW direction if dAMeasure positive
'for Just Basic angle 0 degrees is due East and angle increases clockwise towards South
'dAStart is degrees to start Angle, due East is 0 degrees
'dAMeasure is degrees added (Clockwise) to dAstart for end of arc
rAngleStart = RAD * dAStart rAngleEnd = RAD * dAMeasure + rAngleStart Stepper = RAD / (.1 * arcRadius) 'fixed lastX = xCenter + arcRadius * cos(rAngleStart) lastY = yCenter + arcRadius * sin(rAngleStart) #gr "set ";int(lastX);" ";int(lastY) for rAngle = rAngleStart+Stepper to rAngleEnd step Stepper nextX = xCenter + arcRadius * cos(rAngle) nextY = yCenter + arcRadius * sin(rAngle) #gr "goto ";int(nextX);" ";int(nextY) 'int speeds things up next end sub
These arcs should start same place the Trig functions would graph them and move clockwise as the angle is increased about the center point, same way as the Trig function graphings go!
|
|
|
atan2
Sept 26, 2022 11:45:30 GMT -5
Post by tenochtitlanuk on Sept 26, 2022 11:45:30 GMT -5
Just a nice graphic to confirm I CAN draw arcs reliably! Shades of Sauron's eye, with hgv colouring. This example was created with the arc-ends and radial centre marked. Wihout them and with more arcs drawn you get a smoother but less dramatic pic! This thread also led me to examine the dot.product method. Thanks, friends..
|
|
|
Post by tsh73 on Sept 27, 2022 17:06:45 GMT -5
Ok I'm still at it
To calculate angle between vectors A, B with ATAN2:
get A angle (-pi..pi) calling ATAN2 get B angle (-pi..pi) calling ATAN2 second time subtract. Now, it could end up as (-2pi..2pi) so we should normalize say to 0..2pi (or -pi..pi, whatever)
Now, to get angle by this formula:
calculate vectDotProduct(a$, b$)/ (vectLen(a$)*vectLen(b$)) (that's three user defined functions for me) then take ACS of it !! Alas! ACS returns only angles in 0..pi range Will it be enough? I dunno, might depend on task at hand.
So I used old mine Vector Lib 2d for vector functions run angle A from 0 to 2Pi step pi/6 (variable i, 1st column) run angle B from 0 to 2Pi step pi/6 (variable j, 2nd column) made random vectors with these angles 3rd column - just subtract (and normalise) A, B angles I have 4th - get angle through ATAN2 (happened to == with 3rd, good sign I think) 5th - get angle through ASC and vector functions. (and I had to limit ACS argument because it happen to went out of [-1 .. 1] with some 1e-15 reminder, bringing ACS to error - so it does not run without these checks)
Now, all this is highly inefficient - vectors are stored as strings, I use val(word$())) for getting any numbers Still ACS version works 30% slower. (counted by adding outer loop and removing all else)
'angle between vectors
'using Vector Lib 2d v.01 'vector 2d lib 'vectors as "x y" pairs, to be split by Word$ global pi pi = acs(-1) 'global for atan2, dePi ''
for i = 0 to 11 for j = 0 to 11 a$=vectFromPolar$(1+rnd(0)*5, 2*pi*i/12) b$=vectFromPolar$(1+rnd(0)*5, 2*pi*j/12) da1=(2*pi*j/12)-(2*pi*i/12) 'b$-a$ print using ("##",i);" ";using ("##",j);" "; print dePiNorm$(da1);" "; 'now angle by atan2 a1=vectAngle(a$) b1=vectAngle(b$) da2=b1-a1 'print dePi$(a1);" "; 'so it is from -pi to pi 'print a$;":";b$;","; print dePiNorm$(da2);" "; 'print dePi$(da2);" "; 'now angle by acos 'A = Math.acos ( dot (v1, v2)/ (v1.length ()*v2.length ()) ); dd=vectDotProduct(a$, b$)/ (vectLen(a$)*vectLen(b$)) 'print dd, dd+1;" "; 'somehow it get over -1 if dd<-1 then dd=-1 if dd>1 then dd=1 da3=acs(dd) 'print dePiNorm$(da3) print dePi$(da3)';" "; 'print dd
next next
end '---------------------------------------------
function vect$(x,y) vect$=x;" ";y end function
function vectX(v$) vectX=val(word$(v$,1)) end function
function vectY(v$) vectY=val(word$(v$,2)) end function
function vectLen(v$) x=val(word$(v$,1)) y=val(word$(v$,2)) vectLen=sqr(x*x+y*y) end function
function vectUnit$(v$) x=val(word$(v$,1)) y=val(word$(v$,2)) vectLen=sqr(x*x+y*y) vectUnit$=x/vectLen;" ";y/vectLen end function
function vectAdd$(v1$,v2$) x1=val(word$(v1$,1)) y1=val(word$(v1$,2)) x2=val(word$(v2$,1)) y2=val(word$(v2$,2)) vectAdd$=x1+x2;" ";y1+y2 end function
function vectSub$(v1$,v2$) x1=val(word$(v1$,1)) y1=val(word$(v1$,2)) x2=val(word$(v2$,1)) y2=val(word$(v2$,2)) vectSub$=x1-x2;" ";y1-y2 end function
function vectDotProduct(v1$,v2$) x1=val(word$(v1$,1)) y1=val(word$(v1$,2)) x2=val(word$(v2$,1)) y2=val(word$(v2$,2)) vectDotProduct=x1*x2+y1*y2 end function
function vectScale$(a,v$) 'a * vector v$ x=val(word$(v$,1)) y=val(word$(v$,2)) vectScale$=a*x;" ";a*y end function
function vectTangent$(v$,base$) n$=vectUnit$(base$) vectTangent$=vectScale$(vectDotProduct(n$,v$),n$) end function
function vectNorm$(v$,base$) vectNorm$=vectSub$(v$,vectTangent$(v$,base$)) end function
function vectAngle(v$) x=val(word$(v$,1)) y=val(word$(v$,2)) vectAngle=atan2(y,x) end function
function vectFromPolar$(rho, phi) vectFromPolar$=rho*cos(phi);" ";rho*sin(phi) end function
function vectRotate$(v$,alpha) x=val(word$(v$,1)) y=val(word$(v$,2)) rho=sqr(x*x+y*y) phi=atan2(y,x)+alpha vectRotate$=rho*cos(phi);" ";rho*sin(phi) end function
function dePi$(x) 'pure aestetics 'dePi$=x/pi;"Pi" dePi$=using("##.###",x/pi);"Pi" 'place for sign end function
function normalizeAngle(x) x=x mod (2*pi) if x <0 then x=x+(2*pi) normalizeAngle = x end function
function dePiNorm$(x) 'pure aestetics 'norm z first x=x mod (2*pi) if x <0 then x=x+(2*pi) dePiNorm$=using("##.###",x/pi);"Pi" end function
'--------------------------- function atan2(y,x) 'pi = asn(1) * 2 'global if x <> 0 then arctan = atn(y/x)
select case case x > 0 atan2 = arctan
case y>=0 and x<0 atan2 = pi + arctan
case y<0 and x<0 atan2 = arctan - pi
case y>0 and x=0 atan2 = pi / 2
case y<0 and x=0 atan2 = pi / -2 end select end function
|
|