Binary fractions for cordic?
Graham Stabler
Posts: 2,507
I'm back writing my Cordic object, I have atan, cartesian to polar conversion and asin working, next I do sin, cos, tan, Rsin, Rcos and then on to hyperbolic functions and then perhaps to linear although cordic is not very good for those really. All this is for interger math.
The problem I face is when computing asine, this function mathematically speaking expects a fraction (opposite over hypotenuse) but the general cordic routine the numbers are plugged in to allows for negative numbers. Usually I would consider $FFFF_FFFF to be one with $8000_0000 being a half but instead I have to use only 31 bits for the fraction ignoring the MSB, essentially I have redefined what one is.
When I come to do Sine and Cos then there is the potential for minus fractions so what do I do then, I'm just not sure how to represent numbers for the best. The alternative is to have an argument for one so if you wanted to work out asin(0.5) you could do asin(500,1000) where the 1000 is your new definition of what one actually is, it just can't be $FFFF_FFFF.
But then how do you represent fractions that can be between +/- 1?
I'm just a tad confused at the moment.
Graham
The problem I face is when computing asine, this function mathematically speaking expects a fraction (opposite over hypotenuse) but the general cordic routine the numbers are plugged in to allows for negative numbers. Usually I would consider $FFFF_FFFF to be one with $8000_0000 being a half but instead I have to use only 31 bits for the fraction ignoring the MSB, essentially I have redefined what one is.
When I come to do Sine and Cos then there is the potential for minus fractions so what do I do then, I'm just not sure how to represent numbers for the best. The alternative is to have an argument for one so if you wanted to work out asin(0.5) you could do asin(500,1000) where the 1000 is your new definition of what one actually is, it just can't be $FFFF_FFFF.
But then how do you represent fractions that can be between +/- 1?
I'm just a tad confused at the moment.
Graham
Comments
you could use
top 16 bits is integer part
bottom 16 bits is fractional part.
The top bit is still the sign bit and the fraction is unsigned
int32 represents I.F
with· -32768 <= I <= + 32767
and 0 <= F <= 65535 representing F/65536
Or use the top 2 bits as integer part, and bottom 30 bits as fraction.
to represent values between -1 and +1
-1 = C000_0000
+1 = 4000_0000
0· = 0000_0000
Or use the top bit for sign and other bits for fraction
-1 = FFFF_FFFF (almost -1 but not quite)
0 = 0000_0000
1 = 7FFF_FFFF (almost +1 but not quite)
regards peter
I'll have to get back up to speed about this too, but I'm pretty sure unity should be represented as something like $4000_0000. 30 bits. It has to go + and -, and the code will take advantage of the nice conditional subtract operators in pasm. It can't be $8000_0000, because that would have the sign bit set. And $7FFF_FFFF is very unaesthetic. Also I have found that the accumulation in the successive approximation steps sometimes overflows unity by a small amount, so $4000_0000 is a convenient and aesthetic stand in for unity that leaves room for overflow.
For the angle, I think there too $4000_0000 can stand in for pi/2. I know from sad experience that it is important to leave room for small overflows in the angle accumulator, that is, the angle in successive approximation steps can exceed the value chosen to represent pi/2.
The units for angle do not have to be commensurate or have the same basis as the trig functions of angle. That is $4000_0000 can represent both pi/2 of angle and also unity of the trig functions. (Refer to what Ray Adraka says about the arctan basis function).
You have a point about allowing the user to select the value to represent unity or to allow entries like asin(500,1000). However, I think for achieving the best possible resolution that should be rescaled into a ratio of much larger numbers in order to get the best accuracy from the algorithm. Another problem comes in calculating tan and other functions that are unbounded and need to be represented as either a single number (floating point?) or as a ratio of two integers over a sufficiently large range.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
Tracey, yes $4000_0000 would seem like a neat solution, I've not had problems with overflow because of Chip's idea of prescaling by K rather than post scaling but every little helps.
As far as the angle goes I have just used $FFFF_FFFF to mean 360 and this can't really go wrong. Although cordic does not work over the entire range I can just shift the input argument into the relevant quadrant and then flip the result sign as needed etc. This means you can do a cart-polar(20,135) and get a negative x and positive y value.
I actually use a scaling function to make the numbers big for doing Rsin(theta) etc, I basically shift until the msb is at bit 30 and then shift back at the end.
Tan, hmm good point.
Another thing I'm not so sure about is how useful all of these functions actually are, obviously cartesian to polar transformation is useful and so is the ability to resolve vectors with R.sin(theta) etc but as I'm not really trying to build a calculator will anyone actually find a tan function useful. The log and exp functions will I am sure find use.
Also looking for ideas for a demo program, rather than just plotting some graphs.
Graham
pretty much sure this was the way the floating point co-pro
A la Motorola did it,if you give me 30 mins I'll try and find you
a link!!
Best Regards Ian
:- I'm a dyslexic Satanist,I worship the drivel
Notice that I also mentioned the +1 = 4000_0000 possibility.
In fact, that representation actually has the range [noparse][[/noparse]-2,+2>
so that includes range [noparse][[/noparse]-pi/2,+pi/2].
I need to look in one of the documents from earlier
cordic discussions, but I recall that this representation
was mentioned and it did have some advantages like
representing pi/2 directly.
regards peter
hope it's some use
Nope it's a pdf won't upload
Had this with Old Bit Collector
Had to E-mail him the file
Send me your E-mail and I'll send it you ASAP
"Nope it's a pdf won't upload" - Just change the name extension, or tuck it into a ZIP file.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Beau Schwabe
IC Layout Engineer
Parallax, Inc.
Graham
Using 31 significant digits with the MSB as sign is maybe the simplest approach. Any other would mean using two longs... And would be a pity in that case to not use (some of) those extra bits as significant also.
I've been looking more closely at the floating point methods. FLOATMATH and FLOAT32 unpack the mantissa to put the binary radix to the left of bit 29, just like we were talking about. There is a similar logic to doing that. Calculations can be carried out at higher precision to reduce the propagation of error, and then at the conclusion the mantissa is packed back to 23 bits.
One reason for implementing the functions in CORDIC on the Prop would be to plug the results into floating point. FLOAT32 uses table lookup which has an error bound of 16 bits, so if someone really really needed higher accuracy they could get it via CORDIC.
FLOAT32A has a Tan function that is calculated using a series method, and there are two versions of the function, one for a single argument and one for two arguments as a ratio. Floating point is really good for the unbounded functions like tangent.
You asked about applications for those functions. One that I need is direct calculation of sunset and sunrise times given location on the earth's surface. The equations from the Naval Observatory are heavy spherical trig. I don't recall how much precision it needs, but the equations are pretty stiff. A lot of applications in astronomy are like that. Kinematics? Power factor?
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
One nice thing about the way I am coding this cordic library is that the basic loop is generalized (not my idea), the same one does all circular functions, then there is a different one for all hyperbolic functions, if it seems useul and there is space it could potentially be put on chip which would be excellent for signal processing.
I'm interested in perhaps doing some kinematics myself but most of the stuff I have looked at has required pythagurus rather than trig, delta robots and the like however do need trig. One thing I've realized recently is that if you have an equation with many sins and cos's etc you might as well calculate them at the same time in the same loop, some of the commands at least are common.
Graham
Press the left mouse button and you will see a vector going around and around, its magnitude is defined by the mouse and the angle just sweeps, cordic calculates the x y co-ordinates.
None of the functions are tested fully and there are no linear or hyperbolic functions implemented at all. Don't use the code for anything yet, just thought it might be of interest.
Graham
Any chance for the Atan2 function ?
Scott
Graham
It is not clear to me from a quick glance what you decided about "one". It would help in the heading to have note about the scaling of the input arguments.
Say, atn(arg). What are the units of arg? I think you said you are using $ffff_ffff to represent 2pi, which is the output scaling. From a quick reading of the code, the input ratio is effectively arg/1, where the second argument is simply, #1, no fractional part. So your atn function is computing the arctangent of an integer argument. Am I reading it right?
It is really easier to understand atn2(argy,argx), because the scaling is implied. The angle is a function of the ratio only, and the magnitude should come out in the same units as the arguments. For CORDIC, a wrapper has to scale small input arguments of atn2 up to large values in order to retain numerical precision. (Which your #big subroutine does, entailing no change in computed angle. It does entail rescaling of magnitude after the CORDIC.) The code for ATN2 shoult be pretty much the same as what you have. The #big routine would have to determine which is bigger, x, or y, and adjust to the larger of the two.
Again, I'm in compare and contrast mode with FLOAT32A. There atan is computed as a polynomial in (x-1) / (x+1), for x>0, with seven coefficients stored in a cog ram table. atn2(y,x) is computed by dividing the arguments x:=y/x, and then calling atn(x), and there is an additional wrapper for sign manipulation. The scaling issue is taken care of in floating point by the division. The working variable for the polynomial is, {(x-1) / (x+1), for x>0}, which maps into the interval {-1,+1}. The units of the output angle are commensurate with the units of the input argument. The input argument is of course floating point, "one" being a mantissa with 23 bits of precision. CORDIC could definitely have the advantage of speed, because there are a lot of multiplies and divides in the series solution.
▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
Post Edited (Tracy Allen) : 9/28/2007 5:37:49 PM GMT
At the moment the code sets x to one but this should be the chosen representation of one.
theta = atan(opposite/adjacent), you could supply that ratio in two parts as you say but if adjacent is forced to one then the opposite (y) must equal the supplied ratio un changed. It sort of unwraps the fraction.
Atan2 seems like a no brainer in anycase and is perhaps more sensible anyway.
Graham
Currently $FFFF_FFFF = 0/360 which is OK for some things but I was looking at matlab and it would seem that for atan2 you need a range of +- pi/2 for acos a range of 0 to PI and for a sin also a range of +- pi/2. That makes me wonder if I should have a range of -pi/2 to PI.
But then what about when people come to use these functions, they may want to use them between 0 and 2PI so perhaps the range should be -2PI to 2PI so that would be the top 4 bits for the integer part. This provides a range of -8 to + 8 and still plenty of resolution and conversions to degrees etc will be meaningful.
OK?
Graham
Graham
You already have cartpolar in cordicobject.spin. Isn't that the same as atan2?
Mark
Graham
BTW, I've found your code very helpful in my effort to learn cordic. Thanks for postng it.
Mark