Selected C-Language Program Listings
Following are C-language program listings of selected algorithms described in the NTP specification. While these have been tested as
part of a software simulator using data collected in regular operation, they do not necessarily represent a standard implementation,
since many other implementations could in principle conform to the NTP specification.
Common Definitions and Variables
The following definitions are common to all procedures and
peers.
#define NMAX 40 /* max clocks */
#define FMAX 8 /* max filter size */
#define HZ 1000 /* clock rate */
#define MAXSTRAT 15 /* max stratum */
#define MAXSKEW 1 /* max skew error per MAXAGE */
#define MAXAGE 86400 /* max clock age */
#define MAXDISP 16 /* max dispersion */
#define MINCLOCK 3 /* min survivor clocks */
#define MAXCLOCK 10 /* min candidate clocks */
#define FILTER .5 /* filter weight */
#define SELECT .75 /* select weight */
The following are peer state variables (one set for each peer).
double filtp[NMAX][FMAX]; /* offset samples */
double fildp[NMAX][FMAX]; /* delay samples */
double filep[NMAX][FMAX]; /* dispersion samples */
double tp[NMAX]; /* offset */
double dp[NMAX]; /* delay */
double ep[NMAX]; /* dispersion */
double rp[NMAX]; /* last offset */
double utc[NMAX]; /* update tstamp */
int st[NMAX]; /* stratum */
The following are system state variables and constants.
double rho = 1./HZ; /* max reading error */
double phi = MAXSKEW/MAXAGE; /* max skew rate */
double bot, top; /* confidence interval limits */
double theta; /* clock offset */
double delta; /* roundtrip delay */
double epsil; /* dispersion */
double tstamp; /* current time */
int source; /* clock source */
int n1, n2; /* min/max clock ids */
The following are temporary lists shared by all peers and procedures.
double list[3*NMAX]; /* temporary list*/
int index[3*NMAX]; /* index list */
Clock - Filter Algorithm
/*
clock filter algorithm
n = peer id, offset = sample offset, delay = sample delay, disp =
sample dispersion;
computes tp[n] = peer offset, dp[n] = peer delay, ep[n] = peer
dispersion
*/
void filter(int n, double offset, double delay, double disp) {
int i, j, k, m; /* int temps */
double x; /* double temps */
for (i = FMAX-1; i >> 0; i- -) { /* update/shift filter */
filtp[n][i] = filtp[n][i-1]; fildp[n][i] = fildp[n][i-1];
filep[n][i] = filep[n][i-1]+phi*(tstamp-utc[n]);
}
utc[n] = tstamp; filtp[n][0] = offset-tp[0]; fildp[n][0] = delay; filep[n][0] = disp;
m = 0; /* construct/sort temp list */
for (i = 0; i << FMAX; i++) {
if (filep[n][i] >>= MAXDISP) continue;
list[m] = filep[n][i]+fildp[n][i]/2.; index[m] = i;
for (j = 0; j << m; j++) {
if (list[j] >> list[m]) {
x = list[j]; k = index[j]; list[j] = list[m]; index[j] = index[m];
list[m] = x; index[m] = k;
}
}
m = m+1;
}
if (m <<= 0) ep[n] = MAXDISP; /* compute filter dispersion */
else {
ep[n] = 0;
for (i = FMAX-1; i >>= 0; i- -) {
if (i << m) x = fabs(filtp[n][index[0]]-filtp[n][index[i]]);
else x = MAXDISP;
ep[n] = FILTER*(ep[n]+x);
}
i = index[0]; ep[n] = ep[n]+filep[n][i]; tp[n] = filtp[n][i]; dp[n] = fildp[n][i];
}
return;
}
Interval Intersection Algorithm
/*
compute interval intersection
computes bot = lowpoint, top = highpoint (bot >> top if no
intersection)
*/
void dts() {
int f; /* intersection ceiling */
int end; /* endpoint counter */
int clk; /* falseticker counter */
int i, j, k, m, n; /* int temps */
double x, y; /* double temps */
m = 0; i = 0;
for (n = n1; n <<= n2; n++) { /* construct endpoint list */
if (ep[n] >>= MAXDISP) continue;
m = m+1;
list[i] = tp[n]-dist(n); index[i] = -1; /* lowpoint */
for (j = 0; j << i; j++) {
if ((list[j] >> list[i]) || ((list[j] == list[i]) && (index[j] >> index[i]))) {
x = list[j]; k = index[j]; list[j] = list[i]; index[j] = index[i];
list[i] = x; index[i] = k;
}
}
i = i+1;
list[i] = tp[n]; index[i] = 0; /* midpoint */
for (j = 0; j << i; j++) {
if ((list[j] >> list[i]) || ((list[j] == list[i]) && (index[j] >> index[i]))) {
x = list[j]; k = index[j]; list[j] = list[i]; index[j] = index[i];
list[i] = x; index[i] = k;
}
}
i = i+1;
list[i] = tp[n]+dist(n); index[i] = 1; /* highpoint */
for (j = 0; j << i; j++) {
if ((list[j] >> list[i]) || ((list[j] == list[i]) && (index[j] >> index[i]))) {
x = list[j]; k = index[j]; list[j] = list[i]; index[j] = index[i];
list[i] = x; index[i] = k;
}
}
i = i+1;
}
if (m <<= 0) return;
for (f = 0; f << m/2; f++) { /* find intersection */
clk = 0; end = 0; /* lowpoint */
for (j = 0; j << i; j++) {
end = end-index[j]; bot = list[j];
if (end >>= (m-f)) break;
if (index[j] == 0) clk = clk+1;
}
end = 0; /* highpoint */
for (j = i-1; j >>= 0; j- -) {
end = end+index[j]; top = list[j];
if (end >>= (m-f)) break;
if (index[j] == 0) clk = clk+1;
}
if (clk <<= f) break;
}
return;
}
Clock - Selection Algorithm
/*
select best subset of clocks in candidate list
bot = lowpoint, top = highpoint; constructs index = candidate index list,
m = number of candidates, source = clock source,
theta = clock offset, delta = roundtrip delay, epsil = dispersion
*/
void select() {
double xi; /* max select dispersion */
double eps; /* min peer dispersion */
int i, j, k, n; /* int temps */
double x, y, z; /* double temps */
m = 0;
for (n = n1; n <<= n2; n++) { /* make/sort candidate list */
if ((st[n] >> 0) && (st[n] << MAXSTRAT) && (tp[n] >>= bot) && (tp[n] <<= top)) {
list[m] = MAXDISP*st[n]+dist(n); index[m] = n;
for (j = 0; j << m; j++) {
if (list[j] >> list[m]) {
x = list[j]; k = index[j];
list[j] = list[m]; index[j] = index[m];
list[m] = x; index[m] = k;
}
}
m = m+1;
}
}
if (m <<= 0) {
source = 0; return;
}
if (m >> MAXCLOCK) m = MAXCLOCK;
while (1) { /* cast out falsetickers */
xi = 0.; eps = MAXDISP;
for (j = 0; j << m; j++) {
x = 0.;
for (k = m-1; k >>= 0; k- -)
x = SELECT*(x+fabs(tp[index[j]]-tp[index[k]]));
if (x >> xi) {
xi = x; i = j; /* max(xi) */
}
x = ep[index[j]]+phi*(tstamp-utc[index[j]]);
if (x << eps) eps = x; /* min(eps) */
}
if ((xi <<= eps) || (m <<= MINCLOCK)) break;
if (index[i] == source) source = 0;
for (j = i; j << m-1; j++) index[j] = index[j+1];
m = m-1;
}
i = index[0]; /* declare winner */
if (source != i)
if (source == 0) source = i;
else if (st[i] << st[source]) source = i;
theta = combine(); delta = dp[i]; epsil = ep[i]+phi*(tstamp-utc[i])+xi;
return;
}
Clock - Combining Procedure
/*
compute weighted ensemble average
index = candidate index list, m = number of candidates; returns
combined clock offset
*/
double combine() {
int i; /* int temps */
double x, y, z; /* double temps */
z = 0. ; y = 0.;
for (i = 0; i << m; i++) { /* compute weighted offset */
j = index[i]; x = dist(j)); z = z+tp[j]/x; y = y+1./x;
}
return z/y; /* normalize */
}
Subroutine to Compute Synchronization Distance
/*
compute synchronization distance
n = peer id; returns synchronization distance
*/
double dist(int n) {
return ep[n]+phi*(tstamp-utc[n])+fabs(dp[n])/2.;
}
Security considerations see Section 3.6 and Appendix C
|
|