[ntp:questions] Re: Allan variance plots

David L. Mills mills at udel.edu
Wed Sep 17 18:37:51 UTC 2003


Piotr,

Here is a Matlab program to read loopstats files and construct Allan
deviation plots. You can modify it for another language or file format.
The key is in the first while loop.

Dave

%
% This program computes graphical statistics from
% loopstat files produced by NTP.
%
clear all;
fid = fopen('loopstats', 'r');
[table, count] = fscanf(fid, '%f %f %f %f %f %f %f');
day = table(1:7:count);
sec = table(2:7:count) + (day - day(1)) * 86400;
offset = table(3:7:count);
freq = table(4:7:count);
jitter = table(5:7:count);
wander = table(6:7:count);
poll = table(7:7:count);

y1 = offset; i = 1; d = 16;
while length(y1) >= 10 
	u = diff(y1) / d;
	v = diff(u);
	z1(i) = sqrt(mean(v .* v) / 2);
	m1(i) = d;
	y1 = y1(1:2:length(y1));
	i = i + 1; d = d * 2;
end

loglog(m1, z1 * 1e6, '-k')
axis([1 1e5 1e-4 100]);
xlabel('Time Interval (s)');
ylabel('Allan Deviation (PPM)');
print -dtiff -r600 sim
%
% plot offset
%
clf reset; h = newplot;
set(h, 'FontSize', 16);
set(h, 'Position', [.12 .15 .85 .8]);
plot(sec / 3600, offset * 1e3, '-k');
%axis([0 24 -4 4]);
xlabel('Time (hr)');
ylabel('Offset (ms)');
print -deps offset
print -dtiff -r600 offset
%
xx = abs(offset) * 1e3;
fprintf('mean %.3f median %.3f max %.3f\n', mean(xx), median(xx),
max(xx))
%
% plot jitter
%
clf reset; h = newplot;
set(h, 'FontSize', 16);
set(h, 'Position', [.12 .15 .85 .8]);
plot(sec / 3600, jitter * 1e3, '-k');
%axis([0 24 -4 4]);
xlabel('Time (hr)');
ylabel('Jitter (ms)');
print -deps jitter
print -dtiff -r600 jitter
%
% plot wander
%
clf reset; h = newplot;
set(h, 'FontSize', 16);
set(h, 'Position', [.12 .15 .85 .8]);
plot(sec / 3600, wander, '-k');
%axis([0 24 -4 4]);
xlabel('Time (hr)');
ylabel('Wander (PPM)');
print -deps wander
print -dtiff -r600 wander
%
% plot frequency
%
clf reset; h = newplot;
set(h, 'FontSize', 16);
set(h, 'Position', [.12 .15 .85 .8]);
plot(sec / 3600, freq, '-k');
m = mean(freq);
%axis([0 24 m - .2 m + .2]);
xlabel('Time (hr)');
ylabel('Frequency (PPM)');
print -deps freq
print -dtiff -r600 freq
%
z = sort(abs(offset));
w = 1 - z / max(z);
%
clf reset; h = newplot;
set(h, 'FontSize', 16, 'Position', [.12 .15 .85 .8]);
semilogx(z * 1000, w, '-k');
xlabel('Offset (ms)');
ylabel('P[Offset < y]');
print -deps prob
print -dtiff -r600 prob

Piotr Trojanek wrote:
> 
> Dear NTP experts,
> 
> I have been reading NTP papers for a long time and now I wanted
> to make Allan variance plot of my machine's clock. I have collected
> time differences measured beetwen free running system clock and
> PPS signal from GPS receiver using modified ppsctl.c from PPSkit.
> 
> My data is sampled once per second, so I have offsets like:
> 
> .007540355
> .007513540
> .007486726
> .007460749
> .007433934
> .007412148
> .007380305
> 
> Now I'm confused how to process this data to plots like those in
> NTP papers. I have now about samples from few days of non-stop operation.
> 
> I feed this data directly to dataplot package (NIST software) and use
> 'AV PLOT' command in log-log scale -- result is straight line with
> slope +1, so now I think it is wrong method. I can't figure out the
> right way.
> 
> Thanks in advance for any help!
> 
> --
> Piotr Trojanek



More information about the questions mailing list