Phase continuity of NCO

NCOs can be used to generate frequency sweep signal(chirp). It can sweep smoothly on either side of zero frequency and across zero maintaining continuity of phase. With negative frequency you can reverse the LUT pointer (or invert the frequency through cos/sin switching). However, this feature of phase continuity needs some care.

Once I designed carrier recovery module, based on Costas loop, for a QAM receiver in an FPGA. It worked well in tracking test frequency sweep generated by NCO but it failed at crossing zero frequency. After some work it turned out that NCO output - can range from that nice phase continuous signal to a messy random signal depending on sweep rate relative to its phase resolution.

The following code can be used for testing some cases of NCO sweep. You can set the variable test to 1,2,3,4,5 and see the results, e.g crossing zero smoothly (test 1) or crossing with too much dc(test 2) or phase discontinuity (other tests).

clear all; close all;

test = 1;                              %see comments below
n = 10000;                             %number of test samples

switch test
   case 1, z = linspace(.02,-.01,n);     %gentle zero crossing 
   case 2, z = linspace(.0005,-.0005,n); %excess dc at zero crossing
   case 3, z = randsrc(1,n,[.25 -.25]);  %random sweep, creates messy signal
   case 4, z = randsrc(1,n,[-.25:.001,.25]);  %random sweep, creates messy signal
   case 5, z = ones(1,n/4)*.01; z = [z -z z -z]; %discontinuity at zero crossing
end

%%%%%%%%%%%%%%%%%%%%%%%%% Finite NCO model %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A = 2^15 - 1;                          %max amplitude
M = 2^12;                              %NCO accumulator phase resolution
inc = round(M*z);                      %NCO accumulator phase increment
k = 2^12;                              %lut phase resolution (lut size), 
lsb = log2(M) - log2(k);               %LSBs discarded when addressing
lut = round(A*exp(j*2*pi*(0:k-1)/k));  %lut, one cycle cos/sin data

ptr = 0;
addr = 0;
for i = 1:n
    y(i) = lut(addr+1);                %add 1 for matlab LUT
    ptr = mod(ptr + inc(i), M);        %phase accumulator(ramp)
    addr = round(ptr/2^lsb);           %discard LSBs
    addr(addr >= k) = addr - k;        %check address overflow
end

%display output
plot(real(y),'-');hold
plot(imag(y),'r-');