I'm on a quest, to find an elegant way to convert from rectangular
coordinates to polar coordinates.
What I have, is a pair of 8 bit signed values, for the X and Y coordinates.
I have the radius, which wasn't as ugly as I thought it would be.
Finding the angle in any sort of elegant manner is proving elusive.
Worse, the only manner i know to calculate it is A=csc(X/radius), which
ignores the Y value completely, and has poor resolution where X is near it's
maxima and minima.
The following may help depending on how much accuracy you need.
I have copied this to a few friends who may be able and or willing to
suggest a better alternative.
Have X, Y
Want R & theta (especially theta).
BUT
For small theta X ~ R
and
sin(theta) ~ theta.
But Sin(theta) = Y/R and
Y/R ~ Y/X ~ theta (in radians)
So theta ~ Y/X (!) (in radians)
Lets see if this works in practice.
For say Y/X = 0.1
Theta = 0.09967 radian
Error = +0.33%
= < 1 bit
Not bad
Depends how many bits you want.
A quick Excel simulation gives over 8 bits accuracy up to 6.3
degrees, over 7 bits at 8.5 degrees, over 6 bits at 12.4 degrees,
(Down to 2.2 bits at 45 degrees :-)).
This works fairly well for low or high X/Y (swap X and Y in the
latter case)
If other methods produce better results for higher angles you could
perhaps mix the two depending on the angle.
The error in the basic method increases in a smooth manner and may be
able to be compensated for empirically depending on what mathematical
functions you are able or willing to implement. Basic multiplication
would probably allow
Using only division and multiplication, a bit of empirical playing
gives
Theta predicted = X/Y -0.2(X/Y)^2
which gives an accuracy of -3.45% at 18 degrees (improves to be zero
at 40 degrees) and then increases positively to +1.8% at 45 degrees.
Past 45 degrees use 90-theta.
Below 17 degrees the simple estimate of theta = Y/X is more accurate.
regards
Russell McMahon
PS - rough EXCEL spreadsheet available if required.
(Won't post to list - 100K .XLS, 20K .ZIP)
Need some notion of optimization/resolution/speed criteria.
Interpolation, <abbreviated>series expansion and CORDICS
seem to be the obvious candidates. Maybe we should 'choose sides'
and publish all, targeted toward the 8bit territory. I'll do
(at least) one of them. _I guess I mean the elementary trig
functions to solve just this range of problems_.
> I'm on a quest, to find an elegant way to convert from rectangular
> coordinates to polar coordinates.
>
> What I have, is a pair of 8 bit signed values, for the X and Y coordinates.
> I have the radius, which wasn't as ugly as I thought it would be.
> Finding the angle in any sort of elegant manner is proving elusive.
> Worse, the only manner i know to calculate it is A=csc(X/radius), which
> ignores the Y value completely, and has poor resolution where X is near it's
> maxima and minima.
>
> Anybody know of a way to crack this nut?
> TIA!
> Need some notion of optimization/resolution/speed criteria.
> Interpolation, <abbreviated>series expansion and CORDICS
> seem to be the obvious candidates. Maybe we should 'choose sides'
> and publish all, targeted toward the 8bit territory. I'll do
> (at least) one of them. _I guess I mean the elementary trig
> functions to solve just this range of problems_.
The problem I'm trying to solve, is to "de-noise" an input that gives me an
8 bit "bearing" value. (0-255) over a full circle of rotation.
I may only get a very short burst of samples, (16 is my minimum) or a lot
more, and I intend to use all that I can get.
Each sample is converted to rectangular notation (by sin/cos table lookup)
and I average in four layers of 16 sample pairs, down to a final XY pair.
In the first layer, I use a radius of "1", scaled to 127 to keep the X and Y
values within signed 8 bit values. When the first layer bearings are
averaged into a single entry in the second layer, (also true for 2>3 and
3>4) the X and Y values then carry radus information that will end up as
<=1. If the bearing is very noisy, then the radius will be rather a lot less
than 1, and if quiet, it will closely approach 1.
The final pair are signed bytes, which represent sin and cos on a circle of
radius which is almost certainly not 1. The conventional operation of
A=csc(X/Radius) performs the scaling of X out to the unit circle.
Speed, in my application isn't a huge issue, all the averaging is done
between samples at 300-4000 samples/sec, but this final step is done when
either the input signal stops, or when I have 65000 samples collected, so I
can afford to miss a few while I process.
>> For small theta X ~ R
> One thing, what do you mean by ~ here?
"X is approximately equal to R".
I suspect that table lookup (+/- interpolation) will be the easiest
way to get the first octant arctangent, though Russell's Taylor series
with some optimization will probably perform the same task.
His comments suggest to me that a sparse table lookup and
interpolation may be the most concise. For small Y, you will indeed be
using theta=Y/X, with slightly different slopes above 15 degrees (pi/12
radians). How well would four slope segments approximate?
The radius remains irrelevant to determination of angle, the only
requirement being to scale Y and X at any step to avoid overflow. In
the process, an "binary order of magnitude" indication of radius may be
tallied as a quality indicator.
Russell was muddling X and Y in his explanation. It doesn't matter
of course as what we mean is always "the ratio of Y to X or vice-versa,
being that which is between 0 and 1" and subsequently corrected by
addition to or subtraction from a multiple of pi/2 to locate it in the
appropriate octant.
--
Cheers,
Paul B.
Paul B. Webster VK2BZC wrote:
>
> Dave VanHorn wrote:
>
> >> For small theta X ~ R
> > One thing, what do you mean by ~ here?
>
> "X is approximately equal to R".
> <snipped>
> --
> Cheers,
> Paul B.
Funny story I've heard.
Apparently, one young engineer once decided to REALLY simplify sin(x),
when writing this function in assembly, they decided to just simply
return, not the tiny small x which is approximately the correct answer,
but, instead, 0 (Well, we DID say x is small, right? so it's close to
0, right?)
The craft this sin(x) approximation was supposed to be controlling,
had bad problems with the rudder never deflecting, for "some strange
mysterious reason". <VBG>
Be careful not to OVER-optimize, I guess is the moral of this
story... Software testing can be a good idea, huh?
> Have X, Y
> Want R & theta (especially theta).
> BUT
> For small theta X ~ R
> and
> sin(theta) ~ theta.
> But Sin(theta) = Y/R and
> Y/R ~ Y/X ~ theta (in radians)
> So theta ~ Y/X (!) (in radians)
> Lets see if this works in practice.
I keep tripping over this point, that X~R.
In this application, R isn't fixed, it could be any value between nearly
zero (I wouldn't bother with the rest of the calculation in that case) to
nearly 1
Drastically simplified case: (Normally a minimum of 16 angles would be
averaged)
If I have two input angles, say 5 and 175 degrees, then the averaged angle
would be 90 degrees (amazing how hard the math, and how easy the picture)
but the radius would be very short. The classical approach divides X/R,
and so it performs the scaling, but then I'm required to find the Cosecant
(or some books say Atan)
Before anyone runs off down a blind alley on this, also consider the average
of say 355 and 5 degrees, which is NOT 180.
I can get there, but only if I take the radius value and scale the XY back
out to the unit circle, with some loss of accuracy.
BTW: I entered your approach into mathcad, but it keeps bombing with
"internal error" if I change the angles. :-P
You may have misaprehended the problem and are trying to
get the math to solve the _wrong_ thing. I haven't been following this
thread so this may already have been proposed: but the 'average'
of a set of bearing VECTORS is certainly just their vector sum
(normalized to some length) and since you just want bearing=angle
you can ignore the length.
So if you have the bearings in rectangular coordinates just
do the (multibyte) summations and then the ratio-> trig on the two
most significant bytes. (Involves 16 bit division).
It is difficult for me to see how the 'average' of just two opposite
bearings isn't either zero or UNDEF.
I can't think of any way to compute the vector sum of vectors
in polar coordinates that doesn't _really_ convert to rectangular
to do the summation (which is required for the normal math.
definition of 'average')
Your explanation also seems to suggest you want a 'recent'
average i.e. some time averaging as well. You can get this by saving
the two most sig. bytes of the sum which is at least 8-bit accurate
Also if 'noisy' measurements result (somehow) in shorter
vectors then these would contribute proportionately less to
the vector sum and average..
>
> The problem I'm trying to solve, is to "de-noise" an input that gives me an
> 8 bit "bearing" value. (0-255) over a full circle of rotation.
>
> I may only get a very short burst of samples, (16 is my minimum) or a lot
> more, and I intend to use all that I can get.
>
> Each sample is converted to rectangular notation (by sin/cos table lookup)
> and I average in four layers of 16 sample pairs, down to a final XY pair.
> In the first layer, I use a radius of "1", scaled to 127 to keep the X and Y
> values within signed 8 bit values. When the first layer bearings are
> averaged into a single entry in the second layer, (also true for 2>3 and
> 3>4) the X and Y values then carry radus information that will end up as
> <=1. If the bearing is very noisy, then the radius will be rather a lot less
> than 1, and if quiet, it will closely approach 1.
>
> > Have X, Y
> > Want R & theta (especially theta).
> > BUT
> > For small theta X ~ R
> > and
> > sin(theta) ~ theta.
> > But Sin(theta) = Y/R and
> > Y/R ~ Y/X ~ theta (in radians)
> > So theta ~ Y/X (!) (in radians)
> > Lets see if this works in practice.
>
> I keep tripping over this point, that X~R.
> In this application, R isn't fixed, it could be any value between nearly
> zero (I wouldn't bother with the rest of the calculation in that case) to
> nearly 1
>
> Drastically simplified case: (Normally a minimum of 16 angles would be
> averaged)
> If I have two input angles, say 5 and 175 degrees, then the averaged angle
> would be 90 degrees (amazing how hard the math, and how easy the picture)
> but the radius would be very short. The classical approach divides X/R,
> and so it performs the scaling, but then I'm required to find the Cosecant
> (or some books say Atan)
>
> Before anyone runs off down a blind alley on this, also consider the average
> of say 355 and 5 degrees, which is NOT 180.
>
> I can get there, but only if I take the radius value and scale the XY back
> out to the unit circle, with some loss of accuracy.
>
> BTW: I entered your approach into mathcad, but it keeps bombing with
> "internal error" if I change the angles. :-P
If your angle x isn't NEAR 0, then sin(x) won't be near x.
The MacLaurin Series Expansion formula to use for "perfect" accuracy is:
Sin(x) = x - x^3/3! + x^5/5! - x^7/7! + ... (for x in RADIANS here)
For a SMALL x, x^3 and higher exponents are REALLY small thus can be
thrown away <G>, but if you're not talking about angles near 0 radians
(175 degrees doesn't count!), you can NOT simplify this MacLaurin Series
this way - you need to use more terms in the series, 3 or 4 at least.
It sounds to me like you need to expand the series or use a lookup table
or something. A few terms does a decent job. A lookup table's FAST,
though. If you must do FP Math on a PIC:
Hint for series expansion: Calculate X^2 ONCE, then successive terms
are just multiplied by X^2 & divided by some constants.
Hint #2: Pre-calculate the -1/3! and 1/5! etc. constants, then
MULTIPLY, not divide, by these constants <G>
Hint #3: If you can do it in INT's, it is even faster <G>
> You may have misaprehended the problem and are trying to
> get the math to solve the _wrong_ thing.
I could live with that! :)
>I haven't been following this
> thread so this may already have been proposed: but the 'average'
> of a set of bearing VECTORS is certainly just their vector sum
> (normalized to some length) and since you just want bearing=angle
> you can ignore the length.
Not exactly. My incoming data all has a radius of 1, but I need to average
a lot of bearings, and I need an indictor of how much spread there was in
the data. The radius gives me that.
> So if you have the bearings in rectangular coordinates just
> do the (multibyte) summations and then the ratio-> trig on the two
> most significant bytes. (Involves 16 bit division).
That's where I'm at. I have the grand sum as a 16 bit value for X and Y,
made up of
N samples, where N can be >16, and <=FFFF
> It is difficult for me to see how the 'average' of just two opposite
> bearings isn't either zero or UNDEF.
It is, but the radius is zero (or nearly so) and so I toss it.
> Your explanation also seems to suggest you want a 'recent'
> average i.e. some time averaging as well. You can get this by saving
> the two most sig. bytes of the sum which is at least 8-bit accurate
> Also if 'noisy' measurements result (somehow) in shorter
> vectors then these would contribute proportionately less to
> the vector sum and average..
Yup, that's handled.
A new bearing sample gets converted to rect, and stored in the level 1
buffer.
Each entry from this point onward is an XY pair, not necessarily on the unit
circle.
When level 1 becomes full (16 pairs) then I average level 1 and create a
level 2 entry.
....
When level 4 becomes full, (or if I stop getting bearing inputs) I take the
highest level with data in it, and average that, and calculate the radius.
If the radius is small, then the data is "unusable" and is just dumped.
If the radius is large, then I'll need the result bearing, which is what has
me stuck.
I've got time to calculate though, which is good. While I have data coming
in, it may be as slow as 300S/Sec, or as fast as 4000S/Sec, so the
rectangular math is a big winner on speed.
> If your angle x isn't NEAR 0, then sin(x) won't be near x.
> The MacLaurin Series Expansion formula to use for "perfect" accuracy is:
> Sin(x) = x - x^3/3! + x^5/5! - x^7/7! + ... (for x in RADIANS here)
My input rate precludes using anything but table lookup, I've only got about
200uS to do the conversion, plus averaging 16 values, once every 16 samples.
I am using a lookup table to get SIN and COS values for the polar > rect
conversion.
The table values are signed bytes.
FWIW, my angles are in "degree-like units" of 0-255 over a circle.
I'm simply stuck at the final bit of the rectangular to polar conversion,
Angle=csc(X/R). I have X and Y as 16 bit numbers (The grand sum of the Xs
and Ys), and R would be <=127, where 127 is the unit circle. Taking X/R
scales X out to the unit circle anyway, so all (ALL!) I need is a cosecant
function (in my trig book it's ATAN)
I've tried it as a table, but it dosen't work very well, many angles do not
have unique values within the limits of a 256 byte table. I could go bigger,
but that's expensive.
> I'm simply stuck at the final bit of the rectangular to polar
> conversion,
You sure are. I've been *trying* to tell you the answer, and you are
*determined* to ignore it! One more try:
> Angle=csc(X/R). I have X and Y as 16 bit numbers (The grand sum of
> the Xs and Ys), and R would be <=127, where 127 is the unit circle.
Angle=atn(Y/X) or atn(X/Y) depending on where you place 0¡. Note:
** This is the most linear of the available functions**
To evaluate arctangent, you restrict it to an octant by making signs
match and exchanging X and Y to constrain the ratio Y/X or X/Y to be
positive, less than or equal to one. Then you use either a lookup
table, a sparse lookup table with interpolation, or a Taylor series (or
CORDIC). All these use scaled integer maths.
You use table entries or coefficients scaled to your desired result,
in this case 0 to 32 (one eighth of 256). The final angle is determined
by adding or subtracting this from a quadrant (0, 64, 128, 192)
according to the original signs and comparison of X and Y.
> Taking X/R scales X out to the unit circle anyway,
You don't *need* any calculated determination of R. For your purpose,
take the greater of Y and X (as determined as part of the above process)
and that approximates R to the nearest binary order of magnitude. For
a "quality" figure, that's great!
> so all (ALL!) I need is a cosecant function (in my trig book it's
> ATAN)
That's right. Use the book, Luke! It says Arctangent, it *means*
arctangent.
> I've tried it as a table, but it dosen't work very well, many angles
> do not have unique values within the limits of a 256 byte table. I
> could go bigger, but that's expensive.
> You sure are. I've been *trying* to tell you the answer, and you are
> *determined* to ignore it! One more try:
This is hard for me. Apparently not for you. I apologize for my stupidity.
I take on hard projects, so that I learn something doing them. This one has
been rather a stretch in some regards.
> To evaluate arctangent, you restrict it to an octant by making signs
> match and exchanging X and Y to constrain the ratio Y/X or X/Y to be
> positive, less than or equal to one. Then you use either a lookup
> table, a sparse lookup table with interpolation, or a Taylor series (or
> CORDIC). All these use scaled integer maths.
I follow, I think, I will try to code this.
> That's right. Use the book, Luke! It says Arctangent, it *means*
> arctangent.
I get confused when another book uses csc (I assume cosecant) in the same
position, same formula. My formal math education stopped about mid-high
school, and effectively at the 8th grade, (long ugly story, mistakes I made,
teachers that should have been security guards at a donut factory..)
One of us appears to be confused here. The 'radius' of a
bearing has no meaning as far as I can determine. But you must
mean something by it. If I think of this as a navigation problem
(I used to be an offshore sailor, presently into robot nav.) so
I am used to plotting 'bearing' and distance (= time X rate) to
track position. So you may mean this by 'radius'. And then
'throwing out radius' means ignoring short 'coarses', but
if this is a vector sum that happens (correctly) for free.
If you think of this as a sailboat (or airplane, or something
not constrained to roads) then average bearing (over time) is useful.
But the concept of 'radius of bearing' seems to have little
BEARING <G.> on the problem.
Having gotten SUM(x) and SUM(y) as a (cumulative)
bearing vector, the conversion to angle is conventional.
Dave VanHorn wrote:
> > If your angle x isn't NEAR 0, then sin(x) won't be near x.
Probably me being unclear, I'm pushing my math beyond it's current limits.
>The 'radius' of a
> bearing has no meaning as far as I can determine. But you must
> mean something by it. If I think of this as a navigation problem
> (I used to be an offshore sailor, presently into robot nav.) so
> I am used to plotting 'bearing' and distance (= time X rate) to
> track position. So you may mean this by 'radius'. And then
> 'throwing out radius' means ignoring short 'coarses', but
> if this is a vector sum that happens (correctly) for free.
It's a vector. I have multiple bearings, and as they come in they are just
that.
So I assign them a radius of "1". When you average two vectors, of radius
1, you get a third vector, with a radius of <=1 If the two vectors are
widely separated, then the radius is small, if they are close, the radius is
large, if identical, then the radius is 1.
> If you think of this as a sailboat (or airplane, or something
> not constrained to roads) then average bearing (over time) is useful.
> But the concept of 'radius of bearing' seems to have little
> BEARING <G.> on the problem.
Only as an indication of how confident we are in the bearing.
> Having gotten SUM(x) and SUM(y) as a (cumulative)
> bearing vector, the conversion to angle is conventional.
Once you get past the 'radius' problem then the problem
(as you stated) is to evaluate ATAN (x/y) around the circle.
If you ignore the signs (by proper handling) of x and y then
you need a solution for only one quadrant (90 degrees),
if you further "flip" x and y when |x| > |y| (absolute values)
you end up determining the angle from the vertical OR horizontal
axis (whichever is nearer) i.e solve ATAN (x/y) in an
octant (= 45 degrees). If you want 1 part in 256 around
the unit circle the the 'result' has range of (256/8) = 32
i.e. 5bits accuracy (the other 3 bits tell which octant)
At this point I think I would compute (X/Y) *64, which
since X/Y is in range 0-1, would be in range 0- 63 and
then lookup the result in a 64 value table. So the only
real 'math' is to compute (X/Y)* 64 to produce a <= 6bit
result so takes at most 6 'shift and substract' steps. <Which
is fast even on a Microchip processor <G.>> In particular
NO general purpose mul. or div.
> This is hard for me. Apparently not for you. I apologize for my
> stupidity.
I just get rather frustrated when I'm offering what I have every
confidence is the simple answer, and appear to be ignored.
> I follow, I think, I will try to code this.
The concept of octants and normalization to the first octant rather
than the use of one all-encompassing function was also described by
another poster. That uses simple logic functions instead of maths, and
has side-benefits, such as max(mod(X),mod(Y)) being an adequate
approximation of the Radius as a quality indicator.
> I get confused when another book uses csc (I assume cosecant) in the
> same position, same formula.
Using the *mathematical* convention that 0¡ is the positive X-axis,
Angle=csc(Y/R)
=atn(Y/X)
Swap X and Y if you regard the positive Y-axis as 0¡.
In the first octant (i.e., angle less than 45¡), R is not too
different to X, identical for small values. Insofar as it *is*
different, I realise that the cosecant function is actually more linear
than the arctangent (for this octant - above 45¡ both are hideous).
Since however, this requires intermediate determination of the radius
and *that* requires quite a lot of calculation including a square root,
I am sure it is far easier to use the expression for arctangent
requiring only one operation (division) on the X and Y values
beforehand.
> teachers that should have been security guards at a donut factory..)
Hey, fair go, donut factories need good guards too!
--
Cheers,
Paul B.
It has been pointed out to me that I (and I suspect Dave VanHorn as
well) have confused the terms "cosecant" and "arcsine". My previous
post should properly read:
"Using the *mathematical* convention that 0¡ is the positive X-axis,
Angle=asin(Y/R)
=atan(Y/X)
Swap X and Y if you regard the positive Y-axis as 0¡."
Cosecant is in fact the reciprocal of sine, whilst arcsine, often
denoted sine superscripted -1 or "sine to minus 1" means "the angle
whose sine is".
By way of reminder, all "co" functions equate to the corresponding
function with X and Y reversed and <I>vice versa</I>.
> I'm on a quest, to find an elegant way to convert from rectangular
> coordinates to polar coordinates.
>
> What I have, is a pair of 8 bit signed values, for the X and Y coordinates.
> I have the radius, which wasn't as ugly as I thought it would be.
> Finding the angle in any sort of elegant manner is proving elusive.
> Worse, the only manner i know to calculate it is A=csc(X/radius), which
> ignores the Y value completely, and has poor resolution where X is near it's
> maxima and minima.
>
As others have commented already, it sounds as though if you want an
arctangent as opposed to a 1/sine (though I think you may've meant
arcsine). As I see it, or I guess saw it, there are two problems you need
to solve: 1) Computing the arctangent and 2) performing the division. I
say saw since that's the way I ended up doing it for Jason Wolfson back in
April. Here's what I wrote (privately) to Jason:
I really don't want to 'spoil your fun'. So if you're not interested in a
solution, then don't look at the code that's below :).
I had a little time this afternoon so I decided to hack my arcsine routine
to create the arctan routine. What you see below does not solve your
problem completely - but it is the major part. This routine will take an 8
bit integer that corresponds to the numerator of a fraction whose
denominator is 256 and find its arctangent. So the input ranges from 0 to
255 which corresponds to 0 to 255/256 = 0.996 . The output for an
arctangent routine that returns a floating point number would be from 0
(atan(0)) to 0.783 (atan(255/256)) radians; or if you prefer, 0 to 44.89
degrees. However, this routine scales the output so that pi/4 radians (or
45 degrees) corresponds to 256. So for the input range of 0 to 255 you get
an output of 0 to 255 ( atan(255/256) * 256 / (pi/4) is about 255). It's
probably a little more interesting to see an intermediate data point or
two:
Anyway, I thought I'd share this in case you're interested. The only thing
that's left is combining the fractional division and the swapping of the x
and y values if y is greater than x (and then subtracting the result from
pi/2 or actually 512 in this case).
Regards,
Scott
;----------------------------------------------------------
;
;arctan (as adapted from the similar arcsin function)
;
; The purpose of this routine is to take the arctan of an
;8-bit number that ranges from 0 < x < 255/256. In other
;words, the input, x, is an unsigned fraction whose implicit
;divisor is 256.
; The output is in a conveniently awkward format of binary
;radians (brads?). The output corresponds to the range of zero
;to pi/4 for the normal arctan function. Specifically, this
;algorithm computes:
;
; arctan(x) = real_arctan(x/256) * 256 / (pi/4)
; for 0 <= x <= 255
;
; where, real_arctan returns the real arctan of its argument
;in radians.
;
; The algorithm is a table look-up algorithm plus first order
;linear interpolation. The psuedo code is:
;
;unsigned char arctan(unsigned char x)
;{
; unsigned char i;
;
; i = x >> 4;
; return(arctan[i] + ((arctan[i+1] - arctan[i]) * (x & 0xf))/16);
;}
;
;
arctan
SWAPF x,W
ANDLW 0xf
ADDLW 1
MOVWF temp ;Temporarily store the index
CALL arc_tan_table ;Get a2=atan( (x>>4) + 1)
MOVWF result ;Store temporarily in result
DECF temp,W ;Get the saved index
CALL arc_tan_table ;Get a1=atan( (x>>4) )
SUBWF result,W ;W=a2-a1, This is always positive.
SUBWF result,F ;a1 = a1 - (a1-W) = W
CLRF temp ;Clear the product
CLRC
BTFSC x,0
ADDWF temp,F
RRF temp,F
CLRC
BTFSC x,1
ADDWF temp,F
RRF temp,F
CLRC
BTFSC x,2
ADDWF temp,F
RRF temp,F
CLRC
BTFSC x,3
ADDWF temp,F
RRF temp,W
ADDWF result,F
RETURN
arc_tan_table
ADDWF PCL,F
RETLW 0
RETLW 20 ;atan(1/16) = 3.576deg * 256/45
RETLW 41
RETLW 60
RETLW 80
RETLW 99
RETLW 117
RETLW 134
RETLW 151
RETLW 167
RETLW 182
RETLW 196
RETLW 210
RETLW 222
RETLW 234
RETLW 245
RETLW 0 ;atan(32/32) = 45deg * 256/45
The other part of the problem is implementing the division
FRAC_DIV:
;-------------------
;Fractional division
;
; Given x,y this routine finds:
; a = 256 * y / x
;
movlw 8 ;number of bits in the result
movwf cnt
clrf a ; the result
movf x,w
L1:
clrc
rlf y,f ;if msb of y is set we know x<y
rlf a,f ;and that the lsb of 'a' should be set
subwf y,f ;But we still need to subtract the
;divisor from the dividend just in
;case y is less than 256.
skpnc ;If y>x, but y<256
bsf a,0 ; we still need to set a:0
btfss a,0 ;If y<x then we shouldn't have
addwf y,f ;done the subtraction
decfsz cnt,f
goto L1
return
It's easy enough to combine these two routines to obtain a 4-quadrant
arctan(y/x) routine. However, you do need to keed in mind that the
arctangent routine posted above is only valid over 1/8th of the unit
circle. To obtain the other 7/8th's you'll need to apply the appropriate
trig identities.
Off the top of my head and in C:
// pic routines:
extern int arctan(int x);
extern int frac_div(int y, int x);
// Untested c code that implements a 4-quadrant arctan function
int arctan(int x, int y)
{
int f, swapped;
int reference_angle;
