% Read raw data and construct Allan variance from rawstats data file
% Data are collected using rawstats filtered to remove all peers but usno_raw.
%
clear all;
max = 100000;
wander = .1;
phase = 1;
expon = random('Exponential', phase, max, 1);
%
% Variance-time
%
x1 = linspace(1, length(expon) - 1, length(expon));
y1 = expon; i = 1; d = 1;
while length(y1) >= 100    
	z1(i) = var(y1);
	m1(i) = d;
    clf reset; h = newplot;
    set(h, 'FontSize', 12);
    set(h, 'Position', [.12 .15 .85 .8]);
    plot(x1 / 1000, y1, '-k')
    axis([0 max/1000 0 18]);
    lag = strcat('Time (ks) Interval=', num2str(d));
    xlabel(lag);
    ylabel('Delay (ms)');
    lag = strcat('expon_', num2str(d));
    print('-dtiff', '-r600', lag)
    x1 = (x1(2:2:length(x1)));
    y1 = (y1(2:2:length(y1)) + y1(1:2:length(y1) - 1)) / 2;
    i = i + 1; d = 2 * d;
end
e = m1 / m1(1);
e = 1 ./ e;
loglog(m1, z1, '-k', m1, e * z1(1))
%axis([1 1e5 1e-4 100]);
xlabel('Time Interval (s)');
ylabel('Variance (s^2)');
print -dtiff -r600 expon_var