networkdelta.c revision 1553
1139749Simp/*-
2578Srgrimes * Copyright (c) 1985, 1993
3578Srgrimes *	The Regents of the University of California.  All rights reserved.
41197Srgrimes *
56604Sache * Redistribution and use in source and binary forms, with or without
61197Srgrimes * modification, are permitted provided that the following conditions
71197Srgrimes * are met:
81197Srgrimes * 1. Redistributions of source code must retain the above copyright
91197Srgrimes *    notice, this list of conditions and the following disclaimer.
10578Srgrimes * 2. Redistributions in binary form must reproduce the above copyright
11578Srgrimes *    notice, this list of conditions and the following disclaimer in the
12578Srgrimes *    documentation and/or other materials provided with the distribution.
13578Srgrimes * 3. All advertising materials mentioning features or use of this software
14578Srgrimes *    must display the following acknowledgement:
15578Srgrimes *	This product includes software developed by the University of
16578Srgrimes *	California, Berkeley and its contributors.
17578Srgrimes * 4. Neither the name of the University nor the names of its contributors
18578Srgrimes *    may be used to endorse or promote products derived from this software
19578Srgrimes *    without specific prior written permission.
20578Srgrimes *
21578Srgrimes * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
22578Srgrimes * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
231197Srgrimes * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
24578Srgrimes * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
25578Srgrimes * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
26578Srgrimes * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
27578Srgrimes * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
28578Srgrimes * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29578Srgrimes * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
30578Srgrimes * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
31578Srgrimes * SUCH DAMAGE.
32578Srgrimes */
33578Srgrimes
34578Srgrimes#ifndef lint
35578Srgrimesstatic char sccsid[] = "@(#)networkdelta.c	8.1 (Berkeley) 6/6/93";
36578Srgrimes#endif /* not lint */
37578Srgrimes
38578Srgrimes#ifdef sgi
39578Srgrimes#ident "$Revision: 1.4 $"
40578Srgrimes#endif
41578Srgrimes
42578Srgrimes#include "globals.h"
43578Srgrimes
44119418Sobrienstatic long median __P((float, float *, long *, long *, unsigned int));
45119418Sobrien
46119418Sobrien/*
47260276Sdim * Compute a corrected date.
48578Srgrimes *	Compute the median of the reasonable differences.  First compute
492056Swollman *	the median of all authorized differences, and then compute the
502056Swollman *	median of all differences that are reasonably close to the first
5161011Speter *	median.
522056Swollman *
5324131Sbde * This differs from the original BSD implementation, which looked for
5460041Sphk *	the largest group of machines with essentially the same date.
552056Swollman *	That assumed that machines with bad clocks would be uniformly
56104545Smdodd *	distributed.  Unfortunately, in real life networks, the distribution
5761011Speter *	of machines is not uniform among models of machines, and the
583816Swollman *	distribution of errors in clocks tends to be quite consistent
59104445Smdodd *	for a given model.  In other words, all model VI Supre Servres
60104445Smdodd *	from GoFast Inc. tend to have about the same error.
61104445Smdodd *	The original BSD implementation would chose the clock of the
627430Sbde *	most common model, and discard all others.
63104445Smdodd *
64104445Smdodd *	Therefore, get best we can do is to try to average over all
65104445Smdodd *	of the machines in the network, while discarding "obviously"
66104445Smdodd *	bad values.
67578Srgrimes */
68115477Sphklong
698375Srgrimesnetworkdelta()
70104445Smdodd{
71104445Smdodd	struct hosttbl *htp;
72104445Smdodd	long med;
738375Srgrimes	long lodelta, hidelta;
748375Srgrimes	long logood, higood;
758375Srgrimes	long x[NHOSTS];
76578Srgrimes	long *xp;
772395Sache	int numdelta;
78578Srgrimes	float eps;
79578Srgrimes
8014654Sache	/*
8114654Sache	 * compute the median of the good values
8214654Sache	 */
8314654Sache	med = 0;
8414654Sache	numdelta = 1;
8514654Sache	xp = &x[0];
8614654Sache	*xp = 0;			/* account for ourself */
8714654Sache	for (htp = self.l_fwd; htp != &self; htp = htp->l_fwd) {
8814654Sache		if (htp->good
89578Srgrimes		    && htp->noanswer == 0
90578Srgrimes		    && htp->delta != HOSTDOWN) {
91578Srgrimes			med += htp->delta;
92578Srgrimes			numdelta++;
93578Srgrimes			*++xp = htp->delta;
94578Srgrimes		}
95578Srgrimes	}
961241Sjkh
971241Sjkh	/*
981241Sjkh	 * If we are the only trusted time keeper, then do not change our
991241Sjkh	 * clock.  There may be another time keeping service active.
1001241Sjkh	 */
1011237Sjkh	if (numdelta == 1)
1022477Sache		return 0;
1032477Sache
1044480Sache	med /= numdelta;
1054480Sache	eps = med - x[0];
1064480Sache	if (trace)
10713874Sache		fprintf(fd, "median of %d values starting at %ld is about ",
10813874Sache			numdelta, med);
10913874Sache	med = median(med, &eps, &x[0], xp+1, VALID_RANGE);
1104480Sache
111578Srgrimes	/*
112578Srgrimes	 * compute the median of all values near the good median
113578Srgrimes	 */
114578Srgrimes	hidelta = med + GOOD_RANGE;
115578Srgrimes	lodelta = med - GOOD_RANGE;
116578Srgrimes	higood = med + VGOOD_RANGE;
117578Srgrimes	logood = med - VGOOD_RANGE;
118578Srgrimes	xp = &x[0];
119104445Smdodd	htp = &self;
12011872Sphk	do {
121104445Smdodd		if (htp->noanswer == 0
12211872Sphk		    && htp->delta >= lodelta
123104445Smdodd		    && htp->delta <= hidelta
124104445Smdodd		    && (htp->good
125104445Smdodd			|| (htp->delta >= logood
126104445Smdodd			    && htp->delta <= higood))) {
127104445Smdodd			*xp++ = htp->delta;
128104445Smdodd		}
129104445Smdodd	} while (&self != (htp = htp->l_fwd));
130104445Smdodd
131104445Smdodd	if (xp == &x[0]) {
132104445Smdodd		if (trace)
133104445Smdodd			fprintf(fd, "nothing close to median %ld\n", med);
134104445Smdodd		return med;
135104445Smdodd	}
136104445Smdodd
137141031Ssobomax	if (xp == &x[1]) {
138141031Ssobomax		if (trace)
139104445Smdodd			fprintf(fd, "only value near median is %ld\n", x[0]);
140104445Smdodd		return x[0];
141104445Smdodd	}
14231016Sphk
143104445Smdodd	if (trace)
14431016Sphk		fprintf(fd, "median of %d values starting at %ld is ",
145104445Smdodd		        xp-&x[0], med);
146104445Smdodd	return median(med, &eps, &x[0], xp, 1);
147104445Smdodd}
148104445Smdodd
149104445Smdodd
150104445Smdodd/*
151104445Smdodd * compute the median of an array of signed integers, using the idea
152104445Smdodd *	in <<Numerical Recipes>>.
153104445Smdodd */
154104445Smdoddstatic long
155104445Smdoddmedian(a, eps_ptr, x, xlim, gnuf)
156130585Sphk	float a;			/* initial guess for the median */
157578Srgrimes	float *eps_ptr;			/* spacing near the median */
158104445Smdodd	long *x, *xlim;			/* the data */
159104445Smdodd	unsigned int gnuf;		/* good enough estimate */
160104445Smdodd{
161104445Smdodd	long *xptr;
162578Srgrimes	float ap = LONG_MAX;		/* bounds on the median */
16337389Sjulian	float am = -LONG_MAX;
164126080Sphk	float aa;
165111815Sphk	int npts;			/* # of points above & below guess */
166111815Sphk	float xp;			/* closet point above the guess */
167111815Sphk	float xm;			/* closet point below the guess */
168111815Sphk	float eps;
169111815Sphk	float dum, sum, sumx;
170111815Sphk	int pass;
171126080Sphk#define AMP	1.5			/* smoothing constants */
17247625Sphk#define AFAC	1.5
17337389Sjulian
174578Srgrimes	eps = *eps_ptr;
175578Srgrimes	if (eps < 1.0) {
176578Srgrimes		eps = -eps;
1776604Sache		if (eps < 1.0)
1786604Sache			eps = 1.0;
1796604Sache	}
1806604Sache
181578Srgrimes	for (pass = 1; ; pass++) {	/* loop over the data */
1822477Sache		sum = 0.0;
1832477Sache		sumx = 0.0;
184578Srgrimes		npts = 0;
185578Srgrimes		xp = LONG_MAX;
1862477Sache		xm = -LONG_MAX;
18710069Sjoerg
188578Srgrimes		for (xptr = x; xptr != xlim; xptr++) {
189104445Smdodd			float xx = *xptr;
190104445Smdodd
191578Srgrimes			dum = xx - a;
192104445Smdodd			if (dum != 0.0) {	/* avoid dividing by 0 */
1938876Srgrimes				if (dum > 0.0) {
194104445Smdodd					npts++;
195104445Smdodd					if (xx < xp)
196106490Smdodd						xp = xx;
197104445Smdodd				} else {
198106490Smdodd					npts--;
199578Srgrimes					if (xx > xm)
200578Srgrimes						xm = xx;
201578Srgrimes					dum = -dum;
202104445Smdodd				}
203578Srgrimes				dum = 1.0/(eps + dum);
2044480Sache				sum += dum;
205104545Smdodd				sumx += xx * dum;
206104445Smdodd			}
207104445Smdodd		}
208104545Smdodd
209104445Smdodd		if (ap-am < gnuf || sum == 0) {
210106490Smdodd			if (trace)
211578Srgrimes				fprintf(fd,
212578Srgrimes			           "%ld in %d passes; early out balance=%d\n",
213104441Smdodd				        (long)a, pass, npts);
214130585Sphk			return a;	/* guess was good enough */
215578Srgrimes		}
216104445Smdodd
217104545Smdodd		aa = (sumx/sum-a)*AMP;
2188876Srgrimes		if (npts >= 2) {	/* guess was too low */
219104445Smdodd			am = a;
2208876Srgrimes			aa = xp + max(0.0, aa);;
221578Srgrimes			if (aa > ap)
222106490Smdodd				aa = (a + ap)/2;
223106490Smdodd
224578Srgrimes		} else if (npts <= -2) {  /* guess was two high */
225578Srgrimes			ap = a;
226106490Smdodd			aa = xm + min(0.0, aa);;
227106490Smdodd			if (aa < am)
228578Srgrimes				aa = (a + am)/2;
229104445Smdodd
230106490Smdodd		} else {
231578Srgrimes			break;		/* got it */
232106490Smdodd		}
233106490Smdodd
2346604Sache		if (a == aa) {
235106490Smdodd			if (trace)
236104445Smdodd				fprintf(fd,
237106490Smdodd				  "%ld in %d passes; force out balance=%d\n",
2386604Sache				        (long)a, pass, npts);
2396604Sache			return a;
2406604Sache		}
2416604Sache		eps = AFAC*abs(aa - a);
242106490Smdodd		*eps_ptr = eps;
243104445Smdodd		a = aa;
244106490Smdodd	}
2452477Sache
246106490Smdodd	if (((x - xlim) % 2) != 0) {    /* even number of points? */
247104445Smdodd		if (npts == 0)		/* yes, return an average */
248106490Smdodd			a = (xp+xm)/2;
2492477Sache		else if (npts > 0)
250106490Smdodd			a =  (a+xp)/2;
251104445Smdodd		else
252106490Smdodd			a = (xm+a)/2;
2536604Sache
254578Srgrimes	} else 	if (npts != 0) {	/* odd number of points */
255111731Sphk		if (npts > 0)
256104445Smdodd			a = xp;
257106490Smdodd		else
258104545Smdodd			a = xm;
259578Srgrimes	}
260106490Smdodd
261106490Smdodd	if (trace)
262106490Smdodd		fprintf(fd, "%ld in %d passes\n", (long)a, pass);
263578Srgrimes	return a;
264104545Smdodd}
265106490Smdodd