I have a floating-point digital signal processing system that operates at a fixed sample rate of $f_s = 32768$ samples per second implemented using an x86-64 processor. Assuming that the DSP system is synchronously locked to whatever matters, what is the best way to implement a digital oscillator at some frequency $f$?
Specifically, I want to generate the signal: $$y(t) = \sin(2\pi f t)$$ where $t=n/f_s$ for sample number $n$.
One idea is to keep track of a vector $(x,y)$ which we rotate by an angle $\Delta\phi = 2\pi f/f_s$ on each clock cycle.
As a Matlab pseudocode implementation (the real implementation is in C):
%% Initialization code f_s = 32768; % sample rate [Hz] f = 19.875; % some constant frequency [Hz] v = [1 0]; % initial condition d_phi = 2*pi * f / f_s; % change in angle per clock cycle % initialize the rotation matrix (only once): R = [cos(d_phi), -sin(d_phi) ; ... sin(d_phi), cos(d_phi)]
Then, on each clock cycle, we rotate the vector around a little bit:
%% in-loop code while (forever), v = R*v; % rotate the vector by d_phi y = v(1); % this is the sine wave we're generating output(y); end
This allows the oscillator to be computed with only 4 multiplications per cycle. However, I'd worry about phase error and amplitude stability. (In simple tests I was surprised that the amplitude did not die or explode immediately--maybe the
sincos instruction guarantees $\sin^2+\cos^2=1$?).
What is the right way to do this?