Let’s z be a 2D point in the space as \(z = x + jy\), if we want to rotate this point a given angle \(\theta\), we get the following expressions:
\[e^{j\theta} \cdot z = \left(\cos{\theta} + j \sin{\theta}\right)\left(x+jy\right) \\ = x\cos{\theta}-y\sin{\theta} + j \left(y \cos{\theta} + x \sin{\theta} \right) \\ = x’ + j y’ \]
Then, for a generic point, the rotation can be expressed as an equation system, where \(x’\) and \(y’\) are the new coordinates, \(\theta\) is the rotation angle and \(x\) and \(y\) are the original coordinates:
\[\begin{bmatrix}
x’\\
y’
\end{bmatrix}=
\begin{bmatrix}
\cos{\theta} & -\sin{\theta}\\
\sin{\theta} & \cos{\theta}
\end{bmatrix}\begin{bmatrix}
x\\
y
\end{bmatrix} \]
This rotation can be coded in MATLAB as:
%% Function to rotate vector
function v = rotate(P, theta)
rot = [cos(theta) -sin(theta);
sin(theta) cos(theta)];
v = rot*P;
end
A possible implementation of the cordic algorithm could be:
%% Clear all previous values
clear all;
%% Define vectors
A = [2;-3];
O = [0;0];
%% Define accuracy
error_limit = 0.01;
%% Initialize variables
start_angle = pi/2; % 90º
current_angle = start_angle;
acc_angle = 0;
A_rotated = A;
steps = 0;
% Second quadrant (90º - 180º)
if(A_rotated(1) < 0 && A_rotated(2) > 0)
A_rotated = rotate(A_rotated, -pi/2);
acc_angle = pi/2;
% Third quadrant (180º - 270º)
elseif(A_rotated(1) < 0 && A_rotated(2) < 0)
A_rotated = rotate(A_rotated, -pi);
acc_angle = pi;
% Forth quadrant (270º - 360º)
elseif(A_rotated(1) > 0 && A_rotated(2) < 0)
A_rotated = rotate(A_rotated, -3*pi/2);
acc_angle = 3*pi/2;
end
%% Compute angle
% Keep rotating while error is too high
while(abs(A_rotated(2)) > error_limit)
% Represent current vector
quiver(0, 0,A_rotated(1), A_rotated(2));
% Keep previous vectors
hold on;
% Decrease angle rotation
current_angle = current_angle/2;
% (For debugging purposes)
current_angle_deg = current_angle*180/pi;
% Save current error
error = A_rotated(2);
% If y coordinate is still positive
if(error > 0)
% Rotate again conterclockwise
A_rotated = rotate(A_rotated, -1*current_angle);
% Accumulate rotated angle
acc_angle = acc_angle + current_angle;
% If y coordinate is negative
else
% Rotate vector clockwise
A_rotated = rotate(A_rotated, current_angle);
% Substract current angle to the accumulator because we have
% overcome the actual angle value
acc_angle = acc_angle - current_angle;
end
% (For debugging purposes)
acc_angle_deg = acc_angle*180/pi;
% Increase step counter
steps = steps + 1;
end
% Print angle, hypotenuse length and number of steps needed to compute
fprintf('Angle = %f\nHypotenuse = %f\nNumber of steps: %d\n', acc_angle*180/pi, A_rotated(1), steps);
%% Function to rotate vector
function v = rotate(P, theta)
rot = [cos(theta) -sin(theta);
sin(theta) cos(theta)];
v = rot*P;
end
I have coded an interactive applet to illustrate the algorithm. It has been done using the p5.js library. The error limit has been set to \(0.5\).