fp = fopen('corr.out', 'r');
delta = fscanf(fp, '%e', 1);
[corr numcorr] = fscanf(fp, '%e');
tau = (0:(numcorr-1)) * delta;

zeroval = corr(1) / delta;
s = sprintf('<J^2>/VkT^2 = %.3e W/m/K/s', zeroval);
fprintf(1, '%s\n', s);

clf;
plot (tau, corr/corr(1));
xlabel ('\tau [s]');
ylabel ('<J(0)J(\tau)>/<J^2>');
title (s);
v = axis;
axis([min(tau) max(tau)/10 0 1]);
CHGY;
CHGX;

CPRINT;
