Parallax Forums

Parallax Forums > Public Forums > BASIC Stamp > Tracy Allen cordic examples

Forum Quick Jump

[ << Previous Thread | Next Thread >> ] | Show Newest Post First ]

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/23/2005 6:21 PM (GMT -7)

Hi,

I converted the basic stamp cordic examples at the top of page

http://www.emesystems.com/BS2mathC.htm

to the javelin and get these results

Test program CORDIC

cos( 0.0) = 0.9997   sin( 0.0) = 0.0000
cos( 0.1) = 0.9997   sin( 0.1) = 0.0017
cos( 1.0) = 0.9996   sin( 1.0) = 0.0175
cos( 30.0) = 0.8659   sin( 30.0) = 0.5000
cos( 45.0) = 0.7071   sin( 45.0) = 0.7071
cos( 60.0) = 0.5000   sin( 60.0) = 0.8659
cos( 89.0) = 0.0175   sin( 89.0) = 0.9996
cos( 89.9) = 0.0017   sin( 89.9) = 0.9997
cos( 90.0) = 0.0000   sin( 90.0) = 0.9997
cos( 0.1) = 0.9997   sin( 0.1) = -0.0017
cos( -1.0) = 0.9996   sin( -1.0) = -0.0175
cos(-30.0) = 0.8660   sin(-30.0) = -0.4998
cos(-45.0) = 0.7071   sin(-45.0) = -0.7071
cos(-60.0) = 0.5000   sin(-60.0) = -0.8659
cos(-89.0) = 0.0175   sin(-89.0) = -0.9996
cos(-89.9) = 0.0018   sin(-89.9) = -0.9997
cos(-90.0) = 0.0000   sin(-90.0) = -0.9997
x = 0.0001   y = 0.0001   len = 0.0001   atn = 43.3
x = 0.0010   y = 0.0010   len = 0.0013   atn = 43.3
x = 0.0100   y = 0.0100   len = 0.0141   atn = 45.0
x = 0.1000   y = 0.1000   len = 0.1413   atn = 45.0
x = 1.0000   y = 1.0000   len = 1.4139   atn = 45.0
x = 0.8660   y = 0.5000   len = 0.9997   atn = 30.0
x = -0.8660   y = 0.5000   len = 0.9997   atn = 120.0
x = -0.8660   y = -0.5000   len = 0.9997   atn = -120.0
x = 0.8660   y = -0.5000   len = 0.9997   atn = -30.0
x = 0.5000   y = 0.8660   len = 0.9997   atn = 60.0
x = -0.5000   y = 0.8660   len = 0.9997   atn = 150.0
x = -0.5000   y = -0.8660   len = 0.9997   atn = -150.0
x = 0.5000   y = -0.8660   len = 0.9997   atn = -60.0
program finished

All results are ok, except the calculated arctan values +-120.0 and +-150.0
which should be swapped according to the x and y values.

Can anyone confirm this for the basic stamp?

Also, is the calculated arctan value for x=y=0.0001 correct? (should be 45.0)

but as said the error increases for low x and y.
Also, why is the cordic gain (0.607*65536) set to constant 39793 as it should be 39780?

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/24/2005 2:09 AM (GMT -7)

I found the error in cordic-atn.bse

IF quadrant23 THEN                  ' wrapper for 2nd and 3rd quadrants
    IF zs THEN z=z-9000 ELSE z=z+9000
    ENDIF

should be

IF quadrant23 THEN                  ' wrapper for 2nd and 3rd quadrants
    IF zs THEN z=-18000-z ELSE z=18000-z
    ENDIF

Still want to know about the cordic gain and if the accuracy for low x and y
can be improved.

regards peter

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/24/2005 9:23 AM (GMT -7)

Hi Peter,

You've been busy!

Thanks for the heads up about the error for the quadrant23 business. I'll correct that in the program (but not today--too much holiday stuff to do!). One of the first steps in the program is to move points that happen to be in quadrants 2 or 3 over into quadrant 1 and 4 by reflection in the y axis (x=-x). So once the angle is computed it makes perfect sense to reflect the computed angles back through the y axis and that is what your correction does. I don't know how the rotation by 90 degrees got in there. There must have been dust in the afterburner.

As to the Cordic gain, when I first wrote the programs, I took 0.6072 as an approximation, and then when I explained it in the text I rounded it off to 0.607. That accounts for the difference between 39793 and 39780. But neither one is right. A closer approximation is 0.607252935, and for that the Stamp ** factor should be 39797. I'll work that into the program too, but that probably won't make much difference.

I did go into the www page and changed the exposition. Always tinkering.

As to the errors to calculate arctangent when x=y=0.0001. It is due to integer roundoff. The algorithm needs digits to work with internally to carry through the precision. One solution would be to put another wrapper on it, so that if the input numbers are "small", they will be multiplied times a constant at the outset and then rounded down or handled appropriately at the end. It would be kind of ad-hoc floating point, so that the algorithm always works at maximum precision. That way one could represent the angle in that particular case correctly as 45 degrees and also the length of the vector at 0.00007071. There is no way to get the 0.00007071 without scaling it at the outset. Another way to do it is to go to double precision, 32 bits, but call it 24 bit precision and keep 8 throwaway bits in reserve to slosh around internally. On the demo page, I was more interested in just getting it working on the Stamp--without getting into the complications of double precision.

Here is the sequence of values of the product of the {cos(atn(1/2^i)}, to show how fast it converges.
0.707106781186547
0.632455532033676
0.613571991077896
0.608833912517752
0.607648256256168
0.607351770141296
0.607277644093526
0.607259112298893
0.607254479332562
0.607253321089875
0.607253031529134
0.607252959138945
0.607252941041397
0.607252936517010
0.607252935385914
0.607252935103139
0.607252935032446 <-- 16th product on Stamp **39797
0.607252935014772
0.607252935010354
0.607252935009249
0.607252935008973
0.607252935008904
0.607252935008887
0.607252935008883
0.607252935008882
0.607252935008881
0.607252935008881
0.607252935008881
0.607252935008881

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/24/2005 9:46 AM (GMT -7)

Thanks Tracy for the explanantion.
I found the factor 1/1.64676 = 0.607253 in other references, and this yields 39797.
When I use that my calculation for sin(90)=0.9998, sin(45)=0.7071
Increasing to 39798 yields sin(90)=0.9999, sin(45)=0.7072
So by using 39798 I have 0.0001 errors overall, instead of 0.0000 and 0.0002 errors.
Tweaking this constant appears to allow spreeading the errors evenly.
Changing the constant did not change the calculated atn angles.
Here are the results for 39798

Test program CORDIC

cos( 0.0) = 0.9999   sin( 0.0) = 0.0000
cos( 0.1) = 0.9999   sin( 0.1) = 0.0017
cos( 0.4) = 0.9999   sin( 0.4) = 0.0069
cos( 0.5) = 0.9999   sin( 0.5) = 0.0087
cos( 1.0) = 0.9997   sin( 1.0) = 0.0176
cos( 15.0) = 0.9660   sin( 15.0) = 0.2588
cos( 30.0) = 0.8660   sin( 30.0) = 0.5000
cos( 45.0) = 0.7072   sin( 45.0) = 0.7072
cos( 60.0) = 0.5000   sin( 60.0) = 0.8660
cos( 75.0) = 0.2588   sin( 75.0) = 0.9660
cos( 89.0) = 0.0176   sin( 89.0) = 0.9997
cos( 89.4) = 0.0104   sin( 89.4) = 0.9999
cos( 89.5) = 0.0086   sin( 89.5) = 0.9999
cos( 89.9) = 0.0017   sin( 89.9) = 0.9999
cos( 90.0) = 0.0000   sin( 90.0) = 0.9999
cos( 0.1) = 0.9999   sin( 0.1) = -0.0017
cos( 0.4) = 0.9999   sin( 0.4) = -0.0069
cos( 0.5) = 0.9999   sin( 0.5) = -0.0086
cos( -1.0) = 0.9997   sin( -1.0) = -0.0176
cos(-15.0) = 0.9660   sin(-15.0) = -0.2588
cos(-30.0) = 0.8660   sin(-30.0) = -0.5000
cos(-45.0) = 0.7072   sin(-45.0) = -0.7072
cos(-60.0) = 0.5000   sin(-60.0) = -0.8660
cos(-75.0) = 0.2588   sin(-75.0) = -0.9660
cos(-89.0) = 0.0176   sin(-89.0) = -0.9997
cos(-89.4) = 0.0104   sin(-89.4) = -0.9999
cos(-89.5) = 0.0087   sin(-89.5) = -0.9999
cos(-89.9) = 0.0017   sin(-89.9) = -0.9999
cos(-90.0) = 0.0000   sin(-90.0) = -0.9999
x = 0.0001   y = 0.0001   len = 0.0001   atn = 43.3
x = 0.0010   y = 0.0010   len = 0.0013   atn = 43.3
x = 0.0100   y = 0.0100   len = 0.0141   atn = 45.0
x = 0.1000   y = 0.1000   len = 0.1413   atn = 45.0
x = 1.0000   y = 1.0000   len = 1.4141   atn = 45.0
x = 0.8660   y = 0.5000   len = 0.9998   atn = 30.0
x = -0.8660   y = 0.5000   len = 0.9998   atn = 150.0
x = -0.8660   y = -0.5000   len = 0.9998   atn = -150.0
x = 0.8660   y = -0.5000   len = 0.9998   atn = -30.0
x = 0.5000   y = 0.8660   len = 0.9998   atn = 60.0
x = -0.5000   y = 0.8660   len = 0.9998   atn = 120.0
x = -0.5000   y = -0.8660   len = 0.9998   atn = -120.0
x = 0.5000   y = -0.8660   len = 0.9998   atn = -60.0
program finished

I think I will try the extra scaling for small input numbers.

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/24/2005 11:10 AM (GMT -7)

Tracy, I applied the scaling such that while abs(x) and abs(y) are both <=5000,
I double both x and y (10000 is my unit vector length), and at the end I halve
len for each scale step. That really gives great results (error <= 0.0002 in len)
The results
cos( 0.0) = 0.9999   sin( 0.0) = 0.0000
cos( 0.1) = 0.9999   sin( 0.1) = 0.0017
cos( 0.4) = 0.9999   sin( 0.4) = 0.0069
cos( 0.5) = 0.9999   sin( 0.5) = 0.0087
cos( 1.0) = 0.9997   sin( 1.0) = 0.0176
cos( 15.0) = 0.9660   sin( 15.0) = 0.2588
cos( 30.0) = 0.8660   sin( 30.0) = 0.5000
cos( 45.0) = 0.7072   sin( 45.0) = 0.7072
cos( 60.0) = 0.5000   sin( 60.0) = 0.8660
cos( 75.0) = 0.2588   sin( 75.0) = 0.9660
cos( 89.0) = 0.0176   sin( 89.0) = 0.9997
cos( 89.4) = 0.0104   sin( 89.4) = 0.9999
cos( 89.5) = 0.0086   sin( 89.5) = 0.9999
cos( 89.9) = 0.0017   sin( 89.9) = 0.9999
cos( 90.0) = 0.0000   sin( 90.0) = 0.9999
cos(- 0.1) = 0.9999   sin(- 0.1) = -0.0017
cos(- 0.4) = 0.9999   sin(- 0.4) = -0.0069
cos(- 0.5) = 0.9999   sin(- 0.5) = -0.0086
cos(- 1.0) = 0.9997   sin(- 1.0) = -0.0176
cos(-15.0) = 0.9660   sin(-15.0) = -0.2588
cos(-30.0) = 0.8660   sin(-30.0) = -0.5000
cos(-45.0) = 0.7072   sin(-45.0) = -0.7072
cos(-60.0) = 0.5000   sin(-60.0) = -0.8660
cos(-75.0) = 0.2588   sin(-75.0) = -0.9660
cos(-89.0) = 0.0176   sin(-89.0) = -0.9997
cos(-89.4) = 0.0104   sin(-89.4) = -0.9999
cos(-89.5) = 0.0087   sin(-89.5) = -0.9999
cos(-89.9) = 0.0017   sin(-89.9) = -0.9999
cos(-90.0) = 0.0000   sin(-90.0) = -0.9999
x = 0.0001   y = 0.0001   len = 0.0002   atn = 45.0
x = 0.0010   y = 0.0010   len = 0.0015   atn = 45.0
x = 0.0100   y = 0.0100   len = 0.0142   atn = 45.0
x = 0.1000   y = 0.1000   len = 0.1415   atn = 45.0
x = 1.0000   y = 1.0000   len = 1.4141   atn = 45.0
x = 0.8660   y = 0.5000   len = 0.9998   atn = 30.0
x = -0.8660   y = 0.5000   len = 0.9998   atn = 150.0
x = -0.8660   y = -0.5000   len = 0.9998   atn = -150.0
x = 0.8660   y = -0.5000   len = 0.9998   atn = -30.0
x = 0.5000   y = 0.8660   len = 0.9998   atn = 60.0
x = -0.5000   y = 0.8660   len = 0.9998   atn = 120.0
x = -0.5000   y = -0.8660   len = 0.9998   atn = -120.0
x = 0.5000   y = -0.8660   len = 0.9998   atn = -60.0

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/24/2005 12:18 PM (GMT -7)

Tracy, I found that scaling is also necessary for big input values.
If either abs(x) or abs(y) > 19898 than the calculated atn changes
sign or decreases. This makes sense because x and y are scaled by the cordic gain
so 32767/1.64676 = 19898 is highest allowable input value
In my calculation up to 19900 was allowed before showing errors.
So for big input values both x and y must be scaled down.
I have applied that too and the calculated atn is correct, even for big values.

x=100 and y=20000 returns atn =89.7 which is correct.

regards peter

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/24/2005 12:50 PM (GMT -7)

Yeah, great, that makes sense. I'll put that wrapper on the Stamp code too. I'm really glad you're taking an interest in this.

By the way, someone contacted me recently with an updated URL for the Knaust article, and I found a new link to Ingo Cyliax' BASIC Stamp code in the onliine reader for a course at the U. Alberta. Those updates are in the page at
www.emesystems.com/BS2mathC.htm

Have you posted you Javelin code (class?). I know Andy Lindsey was doing something with CORDIC a while ago for the SX or SX/B, too. I'll make sure he sees this.

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/24/2005 1:49 PM (GMT -7)

Tracy, I have not yet posted it.
I attached my test class. I will create a special cordic class with entries
sin(), cos(), tan(), atan2() etc.
I have attached my testclass for you to see how I adapted your code.
My code allows for angle input in 0.1 degree units, atn results
are also in 0.1 degree units (your code specified 0.01 degree units
for atn results but I think that is not necessary)

regards peter

File Attachment :
cordic2.java 6KB (application/octet-stream)

This file has been downloaded 51 time(s).

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/25/2005 7:25 AM (GMT -7)

Tracy, the atan2(y,x) function really behaves weird when x=0
x = 0.0000   y = 0.0000   atan2(y,x) = - 80.1
x = 0.0000   y = 0.0001   atan2(y,x) = - 90.0
x = 0.0000   y = 0.0010   atan2(y,x) =   90.0
x = 0.0000   y = 0.0100   atan2(y,x) =   90.0
x = 0.0000   y = 0.1000   atan2(y,x) = - 90.0
x = 0.0000   y = 1.0000   atan2(y,x) =   90.0
x = 0.0000   y = -0.0001   atan2(y,x) =   90.0
x = 0.0000   y = -0.0010   atan2(y,x) = - 90.0
x = 0.0000   y = -0.0100   atan2(y,x) = - 90.0
x = 0.0000   y = -0.1000   atan2(y,x) = - 90.0
x = 0.0000   y = -1.0000   atan2(y,x) = - 90.0
x = 0.0001   y = 0.0001   atan2(y,x) =   45.0
x = 0.0001   y = 0.0010   atan2(y,x) =   84.3
x = 0.0001   y = 0.0100   atan2(y,x) =   89.4
x = 0.0001   y = 0.1000   atan2(y,x) =   89.9
x = 0.0001   y = 1.0000   atan2(y,x) =   90.0
x = 0.0001   y = -0.0001   atan2(y,x) = - 45.0
x = 0.0001   y = -0.0010   atan2(y,x) = - 84.3
x = 0.0001   y = -0.0100   atan2(y,x) = - 89.4
x = 0.0001   y = -0.1000   atan2(y,x) = - 89.9
x = 0.0001   y = -1.0000   atan2(y,x) = - 90.0
x = -0.0001   y = 0.0001   atan2(y,x) = 135.0
x = -0.0001   y = 0.0010   atan2(y,x) =   95.7
x = -0.0001   y = 0.0100   atan2(y,x) =   90.6
x = -0.0001   y = 0.1000   atan2(y,x) =   90.1
x = -0.0001   y = 1.0000   atan2(y,x) =   90.0
x = -0.0001   y = -0.0001   atan2(y,x) = -135.0
x = -0.0001   y = -0.0010   atan2(y,x) = - 95.7
x = -0.0001   y = -0.0100   atan2(y,x) = - 90.6
x = -0.0001   y = -0.1000   atan2(y,x) = - 90.1
x = -0.0001   y = -1.0000   atan2(y,x) = - 90.0

When x=0 and y=0 it really is a discontinuity and requires special treatment.
When x=0 and y<>0 the result sometimes have correct sign, sometimes incorrect sign.
The absolute value of the result is correct.
So I test for x=0 as a special case and return angles based on sign of y. Then the results are

x = 0.0000   y = 0.0000   atan2(y,x) =   45.0
x = 0.0000   y = 0.0001   atan2(y,x) =   90.0
x = 0.0000   y = 0.0010   atan2(y,x) =   90.0
x = 0.0000   y = 0.0100   atan2(y,x) =   90.0
x = 0.0000   y = 0.1000   atan2(y,x) =   90.0
x = 0.0000   y = 1.0000   atan2(y,x) =   90.0
x = 0.0000   y = -0.0001   atan2(y,x) = - 90.0
x = 0.0000   y = -0.0010   atan2(y,x) = - 90.0
x = 0.0000   y = -0.0100   atan2(y,x) = - 90.0
x = 0.0000   y = -0.1000   atan2(y,x) = - 90.0
x = 0.0000   y = -1.0000   atan2(y,x) = - 90.0
x = 0.0001   y = 0.0001   atan2(y,x) =   45.0
x = 0.0001   y = 0.0010   atan2(y,x) =   84.3
x = 0.0001   y = 0.0100   atan2(y,x) =   89.4
x = 0.0001   y = 0.1000   atan2(y,x) =   89.9
x = 0.0001   y = 1.0000   atan2(y,x) =   90.0
x = 0.0001   y = -0.0001   atan2(y,x) = - 45.0
x = 0.0001   y = -0.0010   atan2(y,x) = - 84.3
x = 0.0001   y = -0.0100   atan2(y,x) = - 89.4
x = 0.0001   y = -0.1000   atan2(y,x) = - 89.9
x = 0.0001   y = -1.0000   atan2(y,x) = - 90.0
x = -0.0001   y = 0.0001   atan2(y,x) = 135.0
x = -0.0001   y = 0.0010   atan2(y,x) =   95.7
x = -0.0001   y = 0.0100   atan2(y,x) =   90.6
x = -0.0001   y = 0.1000   atan2(y,x) =   90.1
x = -0.0001   y = 1.0000   atan2(y,x) =   90.0
x = -0.0001   y = -0.0001   atan2(y,x) = -135.0
x = -0.0001   y = -0.0010   atan2(y,x) = - 95.7
x = -0.0001   y = -0.0100   atan2(y,x) = - 90.6
x = -0.0001   y = -0.1000   atan2(y,x) = - 90.1
x = -0.0001   y = -1.0000   atan2(y,x) = - 90.0

I treat the case x=0 and y=0 just like y=x and so 45.0 degrees is returned.
I will post my class shortly.

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/26/2005 7:17 AM (GMT -7)

Tracy, I have uploaded my class here:
http://groups.yahoo.com/group/JavelinCode/files/Javelin%20Stamp%20IDE/lib/stamp/math/

Final results:
x=cos( 0.0) = 0.9999   y=sin( 0.0) = 0.0000   tan( 0.0) =    0.00   atan2(y,x) =    0.0
x=cos( 0.1) = 0.9999   y=sin( 0.1) = 0.0017   tan( 0.1) =    0.00   atan2(y,x) =    0.1
x=cos( 0.4) = 0.9999   y=sin( 0.4) = 0.0069   tan( 0.4) =    0.01   atan2(y,x) =    0.4
x=cos( 0.5) = 0.9999   y=sin( 0.5) = 0.0087   tan( 0.5) =    0.01   atan2(y,x) =    0.5
x=cos( 1.0) = 0.9997   y=sin( 1.0) = 0.0176   tan( 1.0) =    0.02   atan2(y,x) =    1.0
x=cos( 15.0) = 0.9660   y=sin( 15.0) = 0.2588   tan( 15.0) =    0.27   atan2(y,x) =   15.0
x=cos( 30.0) = 0.8660   y=sin( 30.0) = 0.5000   tan( 30.0) =    0.58   atan2(y,x) =   30.0
x=cos( 45.0) = 0.7072   y=sin( 45.0) = 0.7072   tan( 45.0) =    1.00   atan2(y,x) =   45.0
x=cos( 60.0) = 0.5000   y=sin( 60.0) = 0.8660   tan( 60.0) =    1.73   atan2(y,x) =   60.0
x=cos( 75.0) = 0.2588   y=sin( 75.0) = 0.9660   tan( 75.0) =    3.73   atan2(y,x) =   75.0
x=cos( 89.0) = 0.0176   y=sin( 89.0) = 0.9997   tan( 89.0) =   56.80   atan2(y,x) =   89.0
x=cos( 89.4) = 0.0104   y=sin( 89.4) = 0.9999   tan( 89.4) =   96.14   atan2(y,x) =   89.4
x=cos( 89.5) = 0.0086   y=sin( 89.5) = 0.9999   tan( 89.5) = 114.93   atan2(y,x) =   89.5
x=cos( 89.8) = 0.0034   y=sin( 89.8) = 0.9999   tan( 89.8) = 294.09   atan2(y,x) =   89.8
x=cos( 89.9) = 0.0017   y=sin( 89.9) = 0.9999   tan( 89.9) = 327.67   atan2(y,x) =   89.9
x=cos( 90.0) = 0.0000   y=sin( 90.0) = 0.9999   tan( 90.0) = 327.67   atan2(y,x) =   90.0
x=cos( 90.1) = -0.0017   y=sin( 90.1) = 0.9999   tan( 90.1) = -327.67   atan2(y,x) =   90.1
x=cos( 90.2) = -0.0034   y=sin( 90.2) = 0.9999   tan( 90.2) = -294.08   atan2(y,x) =   90.2
x=cos(- 0.1) = 0.9999   y=sin(- 0.1) = -0.0017   tan(- 0.1) =    0.00   atan2(y,x) = - 0.1
x=cos(- 0.4) = 0.9999   y=sin(- 0.4) = -0.0069   tan(- 0.4) = - 0.01   atan2(y,x) = - 0.4
x=cos(- 0.5) = 0.9999   y=sin(- 0.5) = -0.0086   tan(- 0.5) = - 0.01   atan2(y,x) = - 0.5
x=cos(- 1.0) = 0.9997   y=sin(- 1.0) = -0.0176   tan(- 1.0) = - 0.02   atan2(y,x) = - 1.0
x=cos(-15.0) = 0.9660   y=sin(-15.0) = -0.2588   tan(-15.0) = - 0.27   atan2(y,x) = - 15.0
x=cos(-30.0) = 0.8660   y=sin(-30.0) = -0.5000   tan(-30.0) = - 0.58   atan2(y,x) = - 30.0
x=cos(-45.0) = 0.7072   y=sin(-45.0) = -0.7072   tan(-45.0) = - 1.00   atan2(y,x) = - 45.0
x=cos(-60.0) = 0.5000   y=sin(-60.0) = -0.8660   tan(-60.0) = - 1.73   atan2(y,x) = - 60.0
x=cos(-75.0) = 0.2588   y=sin(-75.0) = -0.9660   tan(-75.0) = - 3.73   atan2(y,x) = - 75.0
x=cos(-89.0) = 0.0176   y=sin(-89.0) = -0.9997   tan(-89.0) = - 57.13   atan2(y,x) = - 89.0
x=cos(-89.4) = 0.0104   y=sin(-89.4) = -0.9999   tan(-89.4) = - 96.14   atan2(y,x) = - 89.4
x=cos(-89.5) = 0.0087   y=sin(-89.5) = -0.9999   tan(-89.5) = -114.93   atan2(y,x) = - 89.5
x=cos(-89.8) = 0.0034   y=sin(-89.8) = -0.9999   tan(-89.8) = -294.08   atan2(y,x) = - 89.8
x=cos(-89.9) = 0.0017   y=sin(-89.9) = -0.9999   tan(-89.9) = -327.67   atan2(y,x) = - 89.9
x=cos(-90.0) = 0.0000   y=sin(-90.0) = -0.9999   tan(-90.0) = -327.67   atan2(y,x) = - 90.0
x=cos(-90.1) = -0.0017   y=sin(-90.1) = -0.9999   tan(-90.1) = 327.67   atan2(y,x) = - 90.1
x=cos(-90.2) = -0.0034   y=sin(-90.2) = -0.9999   tan(-90.2) = 294.08   atan2(y,x) = - 90.2

regards peter

Newzed
Stamp/Stamp Tester

Date Joined Jul 2004
Total Posts : 2716

Posted 12/26/2005 7:23 AM (GMT -7)

Peter, cos(1) = .9998 and sin(1) = .0174. Are these minor differences important?

Sid

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/26/2005 7:34 AM (GMT -7)

Newzed, the absolute error is <= 0.0002 which I find very good considering
all values are calculated and not read from a table. Also note that the atan2
calculates the original angle given the x=cos and y=sin values (which have
an error in the last decimal as you remarked).

regards peter

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/26/2005 12:15 PM (GMT -7)

Hi Sid,
SIN(1) = 0.0174524
COS(1) = 0.9998477

so the rounding is less than 0.0001 in those cases. The algorithm tends to round down.

Peter, I noticed on the Javelin forum that you also have some trig class code there, but I haven't looked at the Yahoo file site. Is that done with standard polynomials? I wonder, how do you find those compare in terms of speed and accuracy with the CORDIC?

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/26/2005 12:41 PM (GMT -7)

Hi Tracy,
The Trig class uses a 90 words table with sine values for whole degrees. In between
whole degrees the sine and cosine are approximated using
lineair interpolation. The accuray is of the same order
as cordic I guess, but I haven't compare the cordic results
with the table results yet.
regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/26/2005 2:06 PM (GMT -7)

Tracy,
Here are the results using the 90 words table
x=cos( 0.0) = 1.0000   y=sin( 0.0) = 0.0000   tan( 0.0) =    0.00   atan2(y,x) =    0.0
x=cos( 0.1) = 0.9999   y=sin( 0.1) = 0.0017   tan( 0.1) =    0.00   atan2(y,x) =    0.1
x=cos( 0.4) = 0.9999   y=sin( 0.4) = 0.0070   tan( 0.4) =    0.01   atan2(y,x) =    0.4
x=cos( 0.5) = 0.9999   y=sin( 0.5) = 0.0087   tan( 0.5) =    0.01   atan2(y,x) =    0.5
x=cos( 1.0) = 0.9998   y=sin( 1.0) = 0.0175   tan( 1.0) =    0.02   atan2(y,x) =    1.0
x=cos( 15.0) = 0.9659   y=sin( 15.0) = 0.2588   tan( 15.0) =    0.27   atan2(y,x) =   15.1
x=cos( 30.0) = 0.8660   y=sin( 30.0) = 0.5000   tan( 30.0) =    0.58   atan2(y,x) =   30.2
x=cos( 45.0) = 0.7071   y=sin( 45.0) = 0.7071   tan( 45.0) =    1.00   atan2(y,x) =   45.0
x=cos( 60.0) = 0.5000   y=sin( 60.0) = 0.8660   tan( 60.0) =    1.73   atan2(y,x) =   59.8
x=cos( 75.0) = 0.2588   y=sin( 75.0) = 0.9659   tan( 75.0) =    3.73   atan2(y,x) =   74.9
x=cos( 89.0) = 0.0175   y=sin( 89.0) = 0.9998   tan( 89.0) =   57.13   atan2(y,x) =   89.0
x=cos( 89.4) = 0.0105   y=sin( 89.4) = 0.9998   tan( 89.4) =   95.22   atan2(y,x) =   89.4
x=cos( 89.5) = 0.0087   y=sin( 89.5) = 0.9999   tan( 89.5) = 114.93   atan2(y,x) =   89.5
x=cos( 89.8) = 0.0035   y=sin( 89.8) = 0.9999   tan( 89.8) = 285.68   atan2(y,x) =   89.8
x=cos( 89.9) = 0.0017   y=sin( 89.9) = 0.9999   tan( 89.9) = 327.67   atan2(y,x) =   89.9
x=cos( 90.0) = 0.0000   y=sin( 90.0) = 1.0000   tan( 90.0) = 327.67   atan2(y,x) =   90.0
x=cos( 90.1) = -0.0017   y=sin( 90.1) = 0.9999   tan( 90.1) = -327.67   atan2(y,x) =   90.1
x=cos( 90.2) = -0.0035   y=sin( 90.2) = 0.9999   tan( 90.2) = -285.68   atan2(y,x) =   90.2
x=cos(- 0.1) = 0.9999   y=sin(- 0.1) = -0.0017   tan(- 0.1) =    0.00   atan2(y,x) = - 0.1
x=cos(- 0.4) = 0.9999   y=sin(- 0.4) = -0.0070   tan(- 0.4) = - 0.01   atan2(y,x) = - 0.4
x=cos(- 0.5) = 0.9999   y=sin(- 0.5) = -0.0087   tan(- 0.5) = - 0.01   atan2(y,x) = - 0.5
x=cos(- 1.0) = 0.9998   y=sin(- 1.0) = -0.0175   tan(- 1.0) = - 0.02   atan2(y,x) = - 1.0
x=cos(-15.0) = 0.9659   y=sin(-15.0) = -0.2588   tan(-15.0) = - 0.27   atan2(y,x) = - 15.1
x=cos(-30.0) = 0.8660   y=sin(-30.0) = -0.5000   tan(-30.0) = - 0.58   atan2(y,x) = - 30.2
x=cos(-45.0) = 0.7071   y=sin(-45.0) = -0.7071   tan(-45.0) = - 1.00   atan2(y,x) = - 45.0
x=cos(-60.0) = 0.5000   y=sin(-60.0) = -0.8660   tan(-60.0) = - 1.73   atan2(y,x) = - 59.8
x=cos(-75.0) = 0.2588   y=sin(-75.0) = -0.9659   tan(-75.0) = - 3.73   atan2(y,x) = - 74.9
x=cos(-89.0) = 0.0175   y=sin(-89.0) = -0.9998   tan(-89.0) = - 57.13   atan2(y,x) = - 89.0
x=cos(-89.4) = 0.0105   y=sin(-89.4) = -0.9998   tan(-89.4) = - 95.22   atan2(y,x) = - 89.4
x=cos(-89.5) = 0.0087   y=sin(-89.5) = -0.9999   tan(-89.5) = -114.93   atan2(y,x) = - 89.5
x=cos(-89.8) = 0.0035   y=sin(-89.8) = -0.9999   tan(-89.8) = -285.68   atan2(y,x) = - 89.8
x=cos(-89.9) = 0.0017   y=sin(-89.9) = -0.9999   tan(-89.9) = -327.67   atan2(y,x) = - 89.9
x=cos(-90.0) = 0.0000   y=sin(-90.0) = -1.0000   tan(-90.0) = -327.67   atan2(y,x) = - 90.0
x=cos(-90.1) = -0.0017   y=sin(-90.1) = -0.9999   tan(-90.1) = 327.67   atan2(y,x) = - 90.1
x=cos(-90.2) = -0.0035   y=sin(-90.2) = -0.9999   tan(-90.2) = 285.68   atan2(y,x) = - 90.2

It shows the sin and cos are slightly more accurate (and thus tan also), but atan2 is slightly worse.
This atan2 is calculated using the approximation

                    (y/x)            1800
    atn = ------------------------ * ----                 0.1 deg units, for |y/x|<1
          1 + 0.280872*(y/x)*(y/x)    pi

So I favor cordic for atan2 calculation (and that also calculates sqrt(x*x+y*y) ).

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/28/2005 7:05 AM (GMT -7)

Tracy, here are some timing results
Trig uses tables and interpolation, cordic uses calculation

Test speed of trig functions
calculating Trig.sin(x) 100 times takes 107 milliseconds
calculating cordic.sin(x) 100 times takes 4662 milliseconds
calculating Trig.atan2(y,x) 100 times takes 4383 milliseconds
calculating cordic.atan2(y,x) 100 times takes 3561 milliseconds
program finished

So the table approach takes far less time for sin and cos, the atan2
takes less time for cordic calculation (on Javelin). I guess the
basic stamp does a bit better due to its builtin ** operator.

regards peter

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/28/2005 11:30 AM (GMT -7)

Hi Peter,

That makes sense. Table lookup, even with interpolation, will always be faster than iterative calculation.

On the other hand, CORDIC has an edge speedwise on polynomials. CORDIC always improves precision on the order of one bit for each iteration, while polynomials (in general) require an increasing number of iterations for each digit of precision.

I'm not sure how SIN, COS, ATN and HYP functions are implemented in the Stamp, (table or CORDIC). They are limited precision, but knowing Chip Gracey, who has a deep appreciation of this stuff, I'd guess they are CORDIC.

Tracy Allen
www.emesystems.com

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/30/2005 1:46 AM (GMT -7)

Hi Peter,

I updated the cordic-atn.bpe program for the BASIC Stamp to correct the bug you found with reflection into quandrant 2 & 3.
www.emesystems.com/BS2mathC.htm
Also
-- Added a wrapper for scaling small values up to the range of 8192 to 16384 (maintaining the ratio y:x) for better accuracy in the calculation, using Stamp NCD operator.
-- Changed the value used to represent 45 degrees from 16384 to 16200. That avoids the numerical overflow error that occurs in the angle accumulator if it passes the 32767 mark. Along with the scaling, that resolves the strange results you found when x=0.
-- Applied the CORDIC gain before going into the iteration. Applying the gain first scales down large numbers so it can work with larger numbers. Now it works okay for numbers up to at least x=23172 and y=23172. The cutoff is that the length of the vector hits 32768 and thus wraps into negative territory in twos complement.

** with regards

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/30/2005 4:36 AM (GMT -7)

Thanks Tracy,

I noticed the loop executes for idx=0 to 15 but there are

17 entries in the atans table. Is that a typo or must the

loop execute for idx=0 to 16?

regards peter

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/30/2005 5:47 AM (GMT -7)

For those interested,
I found this article that deals with approximation of functions.
There are lots of examples and explanations.

regards peter

File Attachment :
numeric analysis.pdf 248KB (application/pdf)

This file has been downloaded 56 time(s).

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/30/2005 10:43 AM (GMT -7)

Right, there need only be 16 entries in the table. I pulled the values out of an excel spreadsheet and copied one extra, which is zero anyway, so it would not affect the result (?) even if the loop extended FOR 0 TO 16.

A lot of the residual error has to do with roundoff errors in the table entries for the binary arctangents. It occurs to me that a little more precision could be had by applying one more wrapper, to map all the angles into within +/- 45 degrees, and then use a larger number to represent 45 degrees. The wrapper would test x and y and exchange them if abs y > abs x, and after the iteration, adjust the angle by reflection in the 45 degree line and then reflection in the y axis, if the respective reflection flags are set. That in lieu of going to double precision.

I had trouble opening the "numeric analysis" file. Neither Acrobat nor Preview (mac) recognized it and the contents page opened as a graphic in bbedit. It looks interesting. Do you have another URL?

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/30/2005 12:34 PM (GMT -7)

Tracy,

The attachement opens perfectly in adobe reader 7.0.5.
Good point about the 45 degrees. The polynome approximation
uses that also in case |y/x|>1.

regards peter

Tracy Allen
Registered Member

Date Joined Jul 2004
Total Posts : 3219

Posted 12/30/2005 12:47 PM (GMT -7)

Okay, it does open fine with Acrobat reader 7. I was trying earlier with Acrobat 5 and with Preview.

It seems to be chapter 3 from
Peter Turner, _Guide to Numerical Analysis_
www.clarkson.edu/mcs/faculty/turner/
Did you find the whole PDF online somewhere?

Tracy Allen
www.emesystems.com

Peter Verkaik
Registered Member

Date Joined Jul 2004
Total Posts : 3985

Posted 12/30/2005 1:10 PM (GMT -7)

I found this after googling. Unfortenately I lost the exact url (it was turner),
but I am sure this was all. I found this the best resource I looked at because
of the eexamples and the explanations about different approaches.
I don't know if the book is still available.

regards peter

[ << Previous Thread | Next Thread >> ]

Forum Information

Currently it is Thursday, July 29, 2010 5:22 PM (GMT -7)
There are a total of 462,441 posts in 62,066 threads.
In the last 3 days there were 90 new threads and 802 reply posts. View Active Threads