Some CORDIC documentation - Page 2 — Parallax Forums

# Some CORDIC documentation

• Posts: 1,181
edited 2006-10-11 20:45
Thanks everyone for the related links.....I'll have to do some reading.....interesting...we learn something new everyday.

James L
• Posts: 13,610
edited 2006-10-11 21:07
Graham Stabler said...
Race you!

Cordic FFT
What a great idea!

In computing an FFT you must multiply your samples by a sine and a cosine. This is the heart of the FFT and it takes 95% of the computing power - the rest is just adds and increments. CORDIC would simplify this central·operation into a single step. Just set X to your sample, clear Y, rotate by the angle, and you've got the two products in X and Y. The CORDIC expansion could be ignored, since it would be systemic.

CORDIC would do wonders for unraveling the·accummulations into angles and magnitudes, too. That would use the steer-Y-to-zero mode, instead of steer-angle-to-zero used to rotate X and Y.

In fact, this brings the FFT down to "understandable" (never mind the butterfly part), since you're dealing more in polar concepts which are usually obscured by such ill-fitting nomenclature as "real" and "imaginary". It took me forever to understand the root concept of·the·Fourier transform·because of how muddled the definition is due to how math is taught/explained. It's no wonder almost nobody can explain WHY it works, since their brain was prematurely clogged with the formality of the definition. Same goes for integrals and derivatives. I bet 99% of the people who teach calculus can't explain how come reducing an exponent by 1 and multiplying by the original exponent produces a derivative (I don't know, myself, but I'm really curious). Obviously, the person who figured it out understood it, with INTENT. Today, almost nobody knows WHY. They just know the recipe. Education ought to be about understanding WHY, not just HOW. It sure would make things easier for me, anyway. Fortunately, we got Paul Baker at Parallax now and he can "interpret" this math stuff pretty well. And Tracy Allen knows nearly everything, no matter how esoteric, or at least how to figure it out.

Here is one outstanding question that not even Tracy has answered (and it might be because there isn't an answer):

····· How can you sum polar coordinates (angle, distance) without converting back and forth·to/from cartesian? Anyone?

Anyway, Graham, this is a great idea!

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.
• Posts: 2,507
edited 2006-10-11 21:15
I don't know why I didn't do this before

Here is the new spreadsheet that will calculate sin(theta) and cos(theta) given theta and on the second sheet will provide Theta and vector length given X and Y

Find attached

Graham
• Posts: 2,507
edited 2006-10-11 21:24
Chip, since someone mentioned FFTs a while back I have been thinking about it, that idea of the multifunction scope, logic analyser and spectrum analyser got me excited. Not to mention all the other ways you can use FFTs, a fastish photodiode and a laser pointer and you could probably measure blood flow (perfusion) by analyzing the spectrum for example.

When I finally got the cordic method I stuck "cordic fft" into google and got even more excited. I've been trying to refresh myself on FFT's, I used to use them a lot but its been a while.

You are dead right about the way things work being lost in the terminology, e has a lot to answer for.

Graham
• Posts: 13,610
edited 2006-10-11 21:47
Jumping the gun here·a bit, I'm thinking about something that is really needed badly:

A SIMPLE METHOD TO SUM POLAR VECTORS

Here's why: in signal processing we are always "slumming it" in cartesian land because of the simplicity of vector summation: you can simply add the X's and add the Y's. Meanwhile, we're having to lookup or generate sine and cosine, and scale each of them. A lot of work! When it comes time to rotate (X,Y), WHOAH! Thank goodness for CORDIC here. If we could reside in polar land, rotation would be just an add to the angle term. The length could be independently manipulated, too. Never any need to deal with a sine or cosine until we were generating the final output. But, because there's no apparent way to sum vectors in the polar domain, we are forced into the mud of cartesian processing.

If there was an efficient algorithm for summing polar vectors, in polar form, signal processing would become very simple to grasp. Signal processing needs polar treatment. Cartesian is an awful fit.·We need to "live" in polar land just to have the clarity of mind to think about signal processing.

I bet there's some undiscovered algorithm that would make this possible.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.

Post Edited (Chip Gracey (Parallax)) : 10/11/2006 10:13:54 PM GMT
• Posts: 13,610
edited 2006-10-12 05:31
Graham Stabler said...
Chip, since someone mentioned FFTs a while back I have been thinking about it, that idea of the multifunction scope, logic analyser and spectrum analyser got me excited. Not to mention all the other ways you can use FFTs, a fastish photodiode and a laser pointer and you could probably measure blood flow (perfusion) by analyzing the spectrum for example.
That's a fascinating idea. How does the spectrum reflect blood flow, though? It seems the beam is wide enough that it would shine through many blood cells which would probably en masse look like band-limited white noise surging with each pump of the heart. Is that how you would tell?

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.
• Posts: 6,562
edited 2006-10-12 07:45
Chip, regarding the polar addition of vectors. There might be a very constrained way to go about it. Here is a diagram that shows two vectors at angles A1 and A2, lengths L1 and L2,. You want to add them to find A3,L3.

The angle labeled a6 depends only on the angles A1 and A2 of the original vectors, not their lengths. However, the angles a7 and a8, and the resultant angle A3, all do depend on the lengths going in. If there is a constraint on the lengths, if for example they are both the same, or they are unit vectors, or they have some fixed ratio, then it is possible to solve it geometrically and the lengths drop out. For example, if L2 = L1, equal lengths, then,

A3 = A2 + (A1 / 2)

{correction: should be A3 = (A2 + A1) / 2}

Angle A3 depends on only A2 and A1, not L2 and L1 (except that the constraint has to be maintained as the vectors rotate). I haven't thought it through for all the quadrants.

The resultant length is still a problem. Observe again from the geometry that angle a6 depends only on A1 and A2, not on the lengths. if the lengths L1 and L2 are equal as above, then

L3 = 2 * L1 * sin (a6 / 2).

I'm sure that is not the general solution you are looking for. But it could be useful for something. I am thinking of the spirograph, where you see the pattern traced out by resultant sum of two vectors, where the ratio of the rotation frequecies and the ratio of the lengths are both constants.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com

Post Edited (Tracy Allen) : 10/12/2006 5:27:10 PM GMT
• Posts: 3,956
edited 2006-10-12 08:17
Tracy,
for equal lengths, shouldn't that be A3 = (A1+A2)/2
after all, adding (1,0) to (0,1) yields (1,1)
(0 + pi/2)/2 = pi/4

regards peter
• Posts: 6,351
edited 2006-10-12 08:18
Amusing Tracy, I've been working on the magnitude part. The magnitude of the addition of two polar vectors is sqrt(ra2 + rb2 + rarbcos(Θab)). The third term turns out to be the dot product of the two vectors.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Paul Baker
Propeller Applications Engineer

Parallax, Inc.
• Posts: 3,956
edited 2006-10-12 08:27
Paul,
shouldn't that be
sqrt(ra2 + rb2 + 2rarbcos(Θab))
because adding 2 vectors of equals length and angle results
in a vector of double length with the same angle.

regards peter
• Posts: 13,610
edited 2006-10-12 08:32
Tracy Allen said...
Chip, regarding the polar addition of vectors. There might be a very constrained way to go about it. Here is a diagram that shows two vectors at angles A1 and A2, lengths L1 and L2,. You want to add them to find A3,L3.

The angle labeled a6 depends only on the angles A1 and A2 of the original vectors, not their lengths. However, the angles a7 and a8, and the resultant angle A3, all do depend on the lengths going in. If there is a constraint on the lengths, if for example they are both the same, or they are unit vectors, or they have some fixed ratio, then it is possible to solve it geometrically and the lengths drop out. For example, if L2 = L1, equal lengths, then,

A3 = A2 + (A1 / 2)

Angle A3 depends on only A2 and A1, not L2 and L1 (except that the constraint has to be maintained as the vectors rotate). I haven't thought it through for all the quadrants.

The resultant length is still a problem. Observe again from the geometry that angle a6 depends only on A1 and A2, not on the lengths. if the lengths L1 and L2 are equal as above, then

L3 = 2 * L1 * sin (a6 / 2).

I'm sure that is not the general solution you are looking for. But it could be useful for something. I am thinking of the spirograph, where you see the pattern traced out by resultant sum of two vectors, where the ratio of the rotation frequecies and the ratio of the lengths are both constants.

Well, maybe it's analogous to the first iteration of some CORDIC-like algorithm which could be an all-case solution. Thanks for thinking about this. I really believe that with enough thinking, we could come up with at least some algorithm, though not necessarily a processor-efficient one. I've searched the web for this solution, but I've seen on many sites "vectors cannot be added in polar form" (period). I pretty clearly remember having some complex formula for this, though, in my high-school physics class where we needed to add force vectors. If we did come up with·some formula we could start to think about how to break it down.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.
• Posts: 6,351
edited 2006-10-12 08:56
It's late, I'll show my math and see if you can figure if I missed something

rc=sqrt((racos(Θa) + rbcos(Θb))2 + (rasin(Θa) + rbsin(Θb))2)

= sqrt(ra2cos2a) + rb2cos2b) + rarbcos(Θa)cos(Θb) + ra2sin2a) + rb2sin2b) + rarbsin(Θa)sin(Θb))

= sqrt(ra2(cos2a) + sin2a)) + rb2(cos2b) + sin2b)) + rarb(cos(Θa)cos(Θb) + sin(Θa)sin(Θb)))

substituting cos2(&#920[noparse];)[/noparse] + sin2(&#920[noparse];)[/noparse] = 1

= sqrt(ra2 + rb2 + rarb(cos(Θa)cos(Θb) + sin(Θa)sin(Θb)))

substituting sin(Θa)sin(Θb) = (cos(Θab) - cos(Θab))/2

and cos(Θa)cos(Θb) = (cos(Θab)·+ cos(Θab))/2

= sqrt(ra2 + rb2 + rarb(cos(Θab)/2 - cos(Θab)/2 + cos(Θab)/2·+ cos(Θab)/2))

and with the final culmination of terms:

= sqrt(ra2 + rb2 + rarbcos(Θab))

Did I miss something somewhere?

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Paul Baker
Propeller Applications Engineer

Parallax, Inc.

Post Edited (Paul Baker (Parallax)) : 10/12/2006 9:02:21 AM GMT
• Posts: 3,956
edited 2006-10-12 09:00
Chip,

assume 2 vectors in the same quadrant
r1.angle(phi1) and r2.angle(phi2) and r1<r2

result = r1.angle(phi1) + r2.angle(phi2)
= r1.angle(phi1) + r1.angle(phi2) + (r2-r1).angle(phi2)
= r1*K1.angle((phi1+phi2)/2) + (r2-r1).angle(phi2)
= r1'.angle(phi1') + r2'.angle(phi2')
for the next iteration reassign r1,phi1,r2 and phi2 and repeat the above
until r2' reaches near zero (this determines the accuracy)

The difficulty is calculating r1'

regards peter

Post Edited (Peter Verkaik) : 10/12/2006 9:25:47 AM GMT
• Posts: 3,956
edited 2006-10-12 09:04
Paul, you missed the factor 2 in working out (a+b)*(a+b) = a*a + 2*a*b + b*b

It's late, I'll show my math and see if you can figure if I missed something

rc=sqrt((racos(Θa) + rbcos(Θb))2 + (rasin(Θa) + rbsin(Θb))2)

= sqrt(ra2cos2a) + rb2cos2b) + 2rarbcos(Θa)cos(Θb) + ra2sin2a) + rb2sin2b) + 2rarbsin(Θa)sin(Θb))

= sqrt(ra2(cos2a) + sin2a)) + rb2(cos2b) + sin2b)) + 2rarb(cos(Θa)cos(Θb) + sin(Θa)sin(Θb))

substituting cos2(&#920[noparse];)[/noparse] + sin2(&#920[noparse];)[/noparse] = 1

= sqrt(ra2 + rb2 + 2rarb(cos(Θa)cos(Θb) + sin(Θa)sin(Θb))

substituting sin(Θa)sin(Θb) = (cos(Θab) - cos(Θab))/2

and cos(Θa)cos(Θb) = (cos(Θab)·+ cos(Θab))/2

= sqrt(ra2 + rb2 + 2rarb(cos(Θab)/2 - cos(Θab)/2 + cos(Θab)/2·+ cos(Θab)/2)

and with the final culmination of terms:

= sqrt(ra2 + rb2 + 2rarbcos(Θab))

regards peter
• Posts: 6,351
edited 2006-10-12 09:08
Ok I see it now, thanks Peter.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Paul Baker
Propeller Applications Engineer

Parallax, Inc.
• Posts: 2,507
edited 2006-10-12 09:18
Chip Gracey (Parallax) said...

That's a fascinating idea. How does the spectrum reflect blood flow, though? It seems the beam is wide enough that it would shine through many blood cells which would probably en masse look like band-limited white noise surging with each pump of the heart. Is that how you would tell?

You get a doppler signal from the interference between coherent light comming from static tissue (unshifted) and moving tissue(the blood). As there is blood traveling at different speeds and in slightly different directions your spectrum is a mess but if you basically find the average frequency then that will increase with blood flow. I'm not sure you see much pulsing in the capilarys of say the finger where they use this sort of technique.

Graham

Post Edited (Graham Stabler) : 10/12/2006 9:23:01 AM GMT
• Posts: 2,507
edited 2006-10-12 09:28
Just realized I uploaded the old version of the spreadsheet the second time, I'll correct that when I get home.

Graham
• Posts: 231
edited 2006-10-12 11:42
Tracy Allen said...
Chip, regarding the polar addition of vectors. There might be a very constrained way to go about it. Here is a diagram that shows two vectors at angles A1 and A2, lengths L1 and L2,.
I had started doodling about the next step--that of relaxing the equal-length constraint.· Can the vectors be scaled to a unit-length, and their angles "adjusted" by a table lookup (or trig·calculations)·to find the intersection of the original vector thru the circle of unit-length radius?· Then just sum the adjusted angles.

But then it started to look like "slumming" again.

Daniel
• Posts: 231
edited 2006-10-12 12:06
Err, I said that wrong from what I was thinking last night. I'm taking the children on a field trip today, but I'll "re-doodle" and see if it makes sense again.

Daniel
• Posts: 231
edited 2006-10-12 14:31
Tracy Allen said...
If there is a constraint on the lengths, if for example they are both the same, or they are unit vectors, or they have some fixed ratio, then it is possible to solve it geometrically and the lengths drop out.
Where my thinking was going·was·that the GCD of the two lengths seem to me to be a scaling factor, leaving the other factors as the fixed ratio.· If one could quickly find the other factor, then applying Tracy's geometric solution for the angle, and Paul's solution for the resultant length.· Apply this iteratively to each additional polar vector.

I might be wrong.

Daniel
• Posts: 2,507
edited 2006-10-12 15:08
I have now pretty much got my algorithm working, I think it had been working for a while actually without me realizing, one problem was my choice of test angle. 90 degrees seemed like a "simple" choice but no, as Xi and Yi are represented as 32 bit fractions they couldn't go negative and so trying to find the solution around 0 or 90 degrees created errors. So I now have an algorithm that will work between 0 and 90 degrees but not at 0 or 90 degrees.

At this point I'm not too sure what to do in terms of the algorithm, as long as I can produce sin/cos in this region then I can produce it for 0-360 degrees by simply making the answer negative or swapping Xi and Yi etc as if you have cos and sin between 0-90 you have everything you need. However if I was to try angles close to 0 or 90 it might also mess up.

As I want to find the x and y components of a vector of length R I'm wondering if ditching the fractional Xi and Yi and going for R*K for the starting vector length might be a better choice. That would allow the variable to be signed.

The other though was to have the algorithm produce the 32bit fraction that represents Sin and Cos and a sign variable, I'd have to keep track of the sign variable myself during the process though.

Thinking about it I think the option of using R*K as the starting length is probably a win-win situation as the only multiplication is a fixed one that has been optimized by chip.

Comments appreciated, I'll be publishing the code so anyone who thinks they might use it should chime in, also wondering how best to present the code, insertable asm "function" or an object.

Graham
• Posts: 2,507
edited 2006-10-12 15:40
OK, so it is now coded up to have an input angle and vector length and it works, however I can now understand comments made about the size of numbers. When using the fraction the numbers only got small at the end of the function but the dX and dY soon get small now unless the starting value of R is large. I'm guessing I could take the input R and shift the bits to make it bigger and then shift the resulting value back again.

More playing to do.

Graham
• Posts: 2,507
edited 2006-10-12 16:50
It works, I just shift the input R until bit 31 goes high, record the number of shifts needed and then shift back at the end.

The input vector still needs to be of a decent size but in my case I can decide on the size.

I'll post code and example later, I'd still like advice on the best way to do the object.

Graham
• Posts: 6,562
edited 2006-10-12 17:24
Peter Verkaik said...
Tracy,
for equal lengths, shouldn't that be A3 = (A1+A2)/2
after all, adding (1,0) to (0,1) yields (1,1)
(0 + pi/2)/2 = pi/4

regards peter

Thanks for catching that, Peter. You must really be awake, while we are burning the midnight oil over here! I like your disproof by counterexample. There is a lot to look at this morning, with the commentary and and Graham's new results. The midnight oil was bubblng.

The reasoning is as follows:

for angle:
a4 = A1 chord through parallel lines
a5 = A2 - a4 = A2 - A1 composite angle
a6 = 180 - a5 = 180 - (A2 - A1) = 180 + A1 - A2 composite angle
a7 + a8 = 180 - a6 = 180 - (180 + A1 - A2) = A2 - A1 sum of angles of triangle
2 * a8 = A2 - A1 isocoles triangle, L1=L2, a7=a8
a8 = (A2 - A1)/2 algebra
A3 = A1 + a8 = A1 + (A2 -A1)/2 = (A1 + A2) / 2
therefore,
A3 = (A1 + A2)/2

for length:
a6 = 180 + A1 - A2
construct perpendicular to L3 bisecting a6 and bisecting L3
then L3/2 = L1 * sin(a6/2) = L1 sin(90 + (A1 - A2)/2)
under the condition that a6<=180 degrees.

That should jive with the formula you are discussing with Paul, but I don't see the connection right away. It is a special case.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com

Post Edited (Tracy Allen) : 10/12/2006 6:15:23 PM GMT
• Posts: 13,610
edited 2006-10-12 17:48
Graham Stabler said...

I think we all did! Paul, Peter, you, and Tracy and I have all got the bug now. I had a hard time concentrating on what I was working on because I was tempted to start working on this (as in, "again",·since I never get too far before getting re-overwhelmed). If we can get something figured out, I would love to build polar vector summing hardware into the next generation of the chip. It would make DSP so simple to program and understand. There will at least be a CORDIC system in the hub which all cogs could use. But a polar vector summer would be the cat's meow (or the bee's knee's).

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.
• Posts: 3,956
edited 2006-10-12 18:57
Chip, Tracy,

I made a diagram to visualize the iteration process I mentioned.
As it turns out, the vectors that need to be added, move towards
each other until they both have the same length and angle.
At that point the resulting length is just the double value, the angle is the resulting
angle.

In the diagram, start with P0 and Q0. I have drawn circles to find the vector
with smallest length to add to smaller vector. The process is then repeated
with the resulting vector and the remaining vector of the larger vector.
I hope·you understand what I mean.
The iteration should stop when both vectors have·the same angle.

regards peter
·

Post Edited (Peter Verkaik) : 10/12/2006 7:50:36 PM GMT
• Posts: 13,610
edited 2006-10-12 19:16
Peter Verkaik said...
Chip, Tracy,

I made a diagram to visualize the iteration process I mentioned.
As it turns out, the vectors that need to be added, move towards
each other until they both have the same length and angle.
At that point the resulting length is just the double value, the angle is the resulting
angle.

In the diagram, start with P0 and Q0. I have drawn circles to find the vector
with smallest length to add to smaller vector. The process is then repeated
with the resulting vector and the remaining vector of the larger vector.
I hope·you understand what I mean.
The iteration should stop when both angles have the same angle.

regards peter
·
Peter,

This seems like a viable concept that would reduce to shifts and adds. This is exciting! I don't understand all the details of what you're doing, but it definitely feels like the right kind of solution. I think you're onto it!!!

Here's Peter's diagram from his last posting:

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔

Chip Gracey
Parallax, Inc.

Post Edited (Chip Gracey (Parallax)) : 10/12/2006 7:20:57 PM GMT
• Posts: 3,956
edited 2006-10-12 19:27
Chip,
When adding 2 vectors of equal length
the result length = r*sqrt(2+ 2cos(Θab))
so you would need to find a cos for each iteration.
I can't tell if that is faster than
convert polar to coordinate
convert coordinate to polar
as I don't have a clue how fast the angle difference reaches near zero.
Also note the diagram is only for 2 vectors in the same quadrant.

regards peter
·
• Posts: 6,562
edited 2006-10-12 19:37
Hi Peter, thanks for the diagram. I'll have to stare at it later, (cuz now I have to focus on my day job). Now it is you burning the midnight oil! You're in Netherlands, right?

A few months ago I had an email exchange about the Given's transform with a person who had seen my write-up and diagrams, and he was investigating the origins of the transforms in the work of the early astronomers. That goes back to the 17th century and all the way back to Ptolemy. Sometimes those early methods have the seed of something new. This person referred me to a diagram he had published on the 3 dimensional transformation equations, noting that the 2 dimensional equations drop right out. Here is that url if anyone is interested...
www.thediagram.com/4_6/moyer.html

Daniel, I'll be interested to see what you come up with. If the ratio of lengths is constant, it should come in as a multiplier. The CORDIC machinery can do straight multiplcaiton and division and square root, in addition to the fancier functions. They all compute in the same amount of time, so I don't see any advantage in getting rid of a trig function, for example, if it is only to be traded for a multiplcation. _Unless_ the multiplier is built into the hardware! And Chip, that would be fantastic to have the CORDIC engine built into the core. Even a "tranditional", CORDIC engine, very professional, awesome.

▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔▔
Tracy Allen
www.emesystems.com
• Posts: 3,956
edited 2006-10-12 20:09
Hi Tracy,
Yes I live in the Netherlands and it is now 22.10 local time, so not quite midnight here.

Let me try to explain the diagram a little better.
We want to add only vectors of equal length.
Because |Q0| < |P0| I used a circle to find the interception (P1) on vector P0
that has the same length as vector Q0.
R1 is the sum vector Q0+P1
The vector P1 must be subtracted from vector P0. I did this by using an identical
circle (with center P0) to find the interception (P2) that is the remaining vector.
So now we have 2 new vectors: R1 and P2 (R0 = Q0 + P0 = R1 + P2)
Repeat the process for vectors R1 and P2.

Initially I thought the remaining vector length would decrease to zero, but the
diagram shows it is the angle difference between the two new vectors that
decreases to zero.

regards peter