Where does this equation come from?

[curious_mind]

I was referring to this article https://ww1.microchip.com/downloads/en/appnotes/01061a.pdf to understand cordic to implement atan function.
After a struggle, concept is just sinking into my brain. At one point I got stuck up. In Page 5, for vectoring mode, the direction di is given as 1 if yi<0 or -1 otherwise. Where does this relation come from? Is it only by observation or any logic behind this?

 

[FvM]

The explanation is in the text on top of the same page.

 

[std_match]

The logic is that the vector is always rotated in the direction that gets yi closer to zero.

 

[curious_mind]

hi

I do not actually get how yi becomes negative in the first quadrant? ( let us say original angle to be estimated is 30 degrees). If I do conventional calculations with cosphi and sinphi, I get all positive values for x and y for CW and CCW rotations

 

[std_match]

In vectoring mode (page 5), if you start with (x, y) values corresponding to 30 degress, yi will be negative after the first rotation, since the first rotation is 45 degrees. Every rotation will bring yi closer to zero.
The only purpose of di is to decide the direction for the next rotation.

 

[curious_mind]

thanks. I tried it and verified through sample c code. It works. I have one last clarification. Why cordic is restricted to -pi/2 to pi/2.

 

[std_match]

thanks. I tried it and verified through sample c code. It works. I have one last clarification. Why cordic is restricted to -pi/2 to pi/2.

That range is enough since you can always use other ranges by simple trigonometric identities.
However, it should be possible to extend the vector algorithm to +/- pi by adding new first rotations arctan(4) and arctan(2)
It gets a little more complicated since you can't get the next rotation direction by only looking at yi.

It should be possible to extend the rotational mode even more, but the efficiency will drop since the rotation angle will not increase much when you go to arctan(8), arctan(16) etc.
It is better to transform the angle into the range -pi/2 to pi/2 and then use CORDIC.

 

[c_mitra]

Where does this relation come from? Is it only by observation or any logic behind this?

It is for good reason they are called circular trigonometric functions. You can also use them for sin and cos functions. But the Pade approximation is simpler and faster.

The standard expansion for atan(x) converges slowly. There are several devices to get the expansion converge faster. The present method is a variation of Newton Raphson method used in the polar coordinates (if I understand correctly).

I am also not aware of this method but this appears simple and elegant.

 

[wwfeldman]

see the Wiki article:

https://en.wikipedia.org/wiki/CORDIC

there are several references at the end

 

[std_match]

It is for good reason they are called circular trigonometric functions. You can also use them for sin and cos functions. But the Pade approximation is simpler and faster.

That is probably only true if you have a hardware multiplier available with enough number of bits. One feature of the CORDIC algorithm is that there is no multiplication involved. Only shift, add, subtract and a very small table are needed.

CORDIC algorithms helped HP when they designed the revolutionary HP35 calculator, which was released in 1972. How HP35 did it's calculations with extremely limited resources is described in detail here:
http://home.citycable.ch/pierrefleur/Jacques-Laporte/HP 35 Saga.htm

 

[c_mitra]

One feature of the CORDIC algorithm is that there is no multiplication involved. Only shift, add, subtract and a very small table are needed.

This is only "apparently" correct because "shifts" are essentially multiplications or divisions by 2. And if you think in terms of binary system, multiplications are nothing but shifts and additions (particularly integer multiplications).

If you recall the days multiplication was taught in high schools: this is not much different from that.

I do not have a ready reference, integer addition can be very fast, integer multiplication is bit slow and integer division is rather slow (all these are done by the CPU and no need for a hardware coprocessor).

And yes, I remember the HP35- it was priced like a computer! But I still have the Sharp calculator from those bygone days: and it still works.

 

[curious_mind]

Can I get sample C code for 16 bit fixed point cordic, which is understandable?

 

[std_match]

There is pseudo C-code in this article:
https://www.edn.com/design/signal-p...gonometric-Functions-using-a-CORDIC-algorithm

 

[curious_mind]

I tried implementing fixed point C code. I am getting wrong answer. What could be the problem. Answer is 0x4279 instead of 0x4000

int atan_cordic_atan_fixed(int x, int y);
int add(long int a , long int b);
int sub(long int a , long int b);
double test,testx,testy;
int tanLUT_fixed[16];
int p,q,pnew,qnew,acc,final_angle;
int count;
int main()
{
/*Cordic atan range -pi/2 to pi/2 corresponds to to -32768 to 32767.*/
//decimal value= (65535 * angle_deg)/180

test=45; //enter angle

tanLUT_fixed[0]=0x4000; //45 degrees
tanLUT_fixed[1]=0x25c8;
tanLUT_fixed[2]=0x13f6;
tanLUT_fixed[3]=0x0a22;
tanLUT_fixed[4]=0x0516;
tanLUT_fixed[5]=0x028b;
tanLUT_fixed[6]=0x0145;
tanLUT_fixed[7]=0x00A2;
tanLUT_fixed[8]=0x0051;
tanLUT_fixed[9]=0x0029;
tanLUT_fixed[10]=0x0014;
tanLUT_fixed[11]=0x000A;
tanLUT_fixed[12]= 0x0005;
tanLUT_fixed[13]= 0x0002;
tanLUT_fixed[14]= 0x0002;
tanLUT_fixed[15]= 0x0001;


testx=cos(.01745*test);
testy=sin(.01745*test);

final_angle= atan_cordic_atan_fixed((int)(testx*32767),(int)(testy*32767));
printf("final angle is %x\n",final_angle);
return(1);
}
int atan_cordic_atan_fixed(int x, int y)
{
unsigned int i;
p=x;
q=y;
acc=0;
pnew=0;
qnew=0;

for(i=0;i<15;i++)
{
if(q>0) // CW
{
pnew=add(p,(q>>i));
qnew=sub(q,(p>>i));
p=pnew;
q=qnew;
acc=add(acc,tanLUT_fixed);
}
else //CCW
{
pnew=sub(p,(q>>i));
qnew=add(q,(p>>i));
p=pnew;
q=qnew;
acc=sub(acc,tanLUT_fixed);
}

}

return (acc);
}

int add(long int a , long int b)
{
long int result;
result= a+b;
if(result>32767)
{
result=32767;
}
else if((result<-32768))
{
result=32768;
}
return(result);
}

int sub(long int a , long int b)
{
long int result;
result= a-b;
if(result>32767)
{
result=32767;
}
else if((result<-32768))
{
result=32768;
}
return(result);
}

 

[FvM]

Check the calculation flow in debugger, pencil & paper or Excel.
You'll notice that p overflows the 16 Bit int range. Use long for p, pnew, q, qnew.

- - - Updated - - -

Or scale input values with An.

 

[curious_mind]

Hi FvM,

No change in the result by changing the datatype

 

[FvM]

You also need to put hands on your add() and sub() functions.

Perhaps you start with a simple thing, test your code with smaller X and Y values and find out when it fails.

 

[c_mitra]

final_angle= atan_cordic_atan_fixed((int)(testx*32767),(int)(te sty*32767));

There is a space present; remove the space and recompile. You cannot return a long int (i.e., result) as a regular int (in the two functions add and sub); the compiler may not grumble but the results will be unpredictable.

 

[curious_mind]

Added the printf in for loop: printf("p,q,accumulator data %d\t%d\t%f\n",(int)p,(int)q,(double)((acc*180.0)/65535));

results are recorded below. The q value does not change sign and hence wrong result. It is strange. Initial values p=23173, q=23166 ( expected value is 45 degrees)

p,q,accumulator data 32767 -7 45.000687
p,q,accumulator data 32767 16376 18.435339
p,q,accumulator data 32767 8185 32.470589
p,q,accumulator data 32767 4090 39.595329
p,q,accumulator data 32767 2043 43.171436
p,q,accumulator data 32767 1020 44.959488
p,q,accumulator data 32767 509 45.852139
p,q,accumulator data 32767 254 46.297092
p,q,accumulator data 32767 127 46.519569
p,q,accumulator data 32767 64 46.632183
p,q,accumulator data 32767 33 46.687115
p,q,accumulator data 32767 18 46.714581
p,q,accumulator data 32767 11 46.728313
p,q,accumulator data 32767 8 46.733807
p,q,accumulator data 32767 7 46.739300

 

[FvM]

These are the intermediate result if you don't clip the values to 16 Bit int



You see why the number range can't be handled in 16 bit