1 /* $NetBSD: networkdelta.c,v 1.13 2021/10/18 14:16:49 christos Exp $ */
2
3 /*-
4 * Copyright (c) 1985, 1993 The Regents of the University of California.
5 * All rights reserved.
6 *
7 * Redistribution and use in source and binary forms, with or without
8 * modification, are permitted provided that the following conditions
9 * are met:
10 * 1. Redistributions of source code must retain the above copyright
11 * notice, this list of conditions and the following disclaimer.
12 * 2. Redistributions in binary form must reproduce the above copyright
13 * notice, this list of conditions and the following disclaimer in the
14 * documentation and/or other materials provided with the distribution.
15 * 3. Neither the name of the University nor the names of its contributors
16 * may be used to endorse or promote products derived from this software
17 * without specific prior written permission.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
20 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
23 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
25 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
26 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
29 * SUCH DAMAGE.
30 */
31
32 #include <sys/cdefs.h>
33 #ifndef lint
34 #if 0
35 static char sccsid[] = "@(#)networkdelta.c 8.3 (Berkeley) 4/27/95";
36 #else
37 __RCSID("$NetBSD: networkdelta.c,v 1.13 2021/10/18 14:16:49 christos Exp $");
38 #endif
39 #endif /* not lint */
40
41 #include "globals.h"
42
43 static long median(float, float*, long*, long*, unsigned int);
44
45 /*
46 * Compute a corrected date.
47 * Compute the median of the reasonable differences. First compute
48 * the median of all authorized differences, and then compute the
49 * median of all differences that are reasonably close to the first
50 * median.
51 *
52 * This differs from the original BSD implementation, which looked for
53 * the largest group of machines with essentially the same date.
54 * That assumed that machines with bad clocks would be uniformly
55 * distributed. Unfortunately, in real life networks, the distribution
56 * of machines is not uniform among models of machines, and the
57 * distribution of errors in clocks tends to be quite consistent
58 * for a given model. In other words, all model VI Supre Servres
59 * from GoFast Inc. tend to have about the same error.
60 * The original BSD implementation would chose the clock of the
61 * most common model, and discard all others.
62 *
63 * Therefore, get best we can do is to try to average over all
64 * of the machines in the network, while discarding "obviously"
65 * bad values.
66 */
67 long
networkdelta(void)68 networkdelta(void)
69 {
70 struct hosttbl *htp;
71 long med;
72 long lodelta, hidelta;
73 long logood, higood;
74 long x[NHOSTS];
75 long *xp;
76 int numdelta;
77 float eps;
78
79 /*
80 * compute the median of the good values
81 */
82 med = 0;
83 numdelta = 1;
84 xp = &x[0];
85 *xp = 0; /* account for ourself */
86 for (htp = self.l_fwd; htp != &self; htp = htp->l_fwd) {
87 if (htp->good
88 && htp->noanswer == 0
89 && htp->delta != HOSTDOWN) {
90 med += htp->delta;
91 numdelta++;
92 *++xp = htp->delta;
93 }
94 }
95
96 /*
97 * If we are the only trusted time keeper, then do not change our
98 * clock. There may be another time keeping service active.
99 */
100 if (numdelta == 1)
101 return 0;
102
103 med /= numdelta;
104 eps = med - x[0];
105 if (trace)
106 fprintf(fd, "median of %d values starting at %ld is about ",
107 numdelta, med);
108 med = median(med, &eps, &x[0], xp+1, VALID_RANGE);
109
110 /*
111 * compute the median of all values near the good median
112 */
113 hidelta = med + GOOD_RANGE;
114 lodelta = med - GOOD_RANGE;
115 higood = med + VGOOD_RANGE;
116 logood = med - VGOOD_RANGE;
117 xp = &x[0];
118 htp = &self;
119 do {
120 if (htp->noanswer == 0
121 && htp->delta >= lodelta
122 && htp->delta <= hidelta
123 && (htp->good
124 || (htp->delta >= logood
125 && htp->delta <= higood))) {
126 *xp++ = htp->delta;
127 }
128 } while (&self != (htp = htp->l_fwd));
129
130 if (xp == &x[0]) {
131 if (trace)
132 fprintf(fd, "nothing close to median %ld\n", med);
133 return med;
134 }
135
136 if (xp == &x[1]) {
137 if (trace)
138 fprintf(fd, "only value near median is %ld\n", x[0]);
139 return x[0];
140 }
141
142 if (trace)
143 fprintf(fd, "median of %ld values starting at %ld is ",
144 (long)(xp - &x[0]), med);
145 return median(med, &eps, &x[0], xp, 1);
146 }
147
148
149 /*
150 * compute the median of an array of signed integers, using the idea
151 * in <<Numerical Recipes>>.
152 */
153 static long
median(float a,float * eps_ptr,long * x,long * xlim,unsigned int gnuf)154 median(float a, /* initial guess for the median */
155 float *eps_ptr, /* spacing near the median */
156 long *x, long *xlim, /* the data */
157 unsigned int gnuf) /* good enough estimate */
158 {
159 long *xptr;
160 float ap = (float)LONG_MAX; /* bounds on the median */
161 float am = (float)-LONG_MAX;
162 float aa;
163 int npts; /* # of points above & below guess */
164 float xp; /* closet point above the guess */
165 float xm; /* closet point below the guess */
166 float eps;
167 float dum, sum, sumx;
168 int pass;
169 #define AMP 1.5 /* smoothing constants */
170 #define AFAC 1.5
171
172 eps = *eps_ptr;
173 if (eps < 1.0) {
174 eps = -eps;
175 if (eps < 1.0)
176 eps = 1.0;
177 }
178
179 for (pass = 1; ; pass++) { /* loop over the data */
180 sum = 0.0;
181 sumx = 0.0;
182 npts = 0;
183 xp = (float)LONG_MAX;
184 xm = (float)-LONG_MAX;
185
186 for (xptr = x; xptr != xlim; xptr++) {
187 float xx = *xptr;
188
189 dum = xx - a;
190 if (dum != 0.0) { /* avoid dividing by 0 */
191 if (dum > 0.0) {
192 npts++;
193 if (xx < xp)
194 xp = xx;
195 } else {
196 npts--;
197 if (xx > xm)
198 xm = xx;
199 dum = -dum;
200 }
201 dum = 1.0/(eps + dum);
202 sum += dum;
203 sumx += xx * dum;
204 }
205 }
206
207 if (ap-am < gnuf || sum == 0) {
208 if (trace)
209 fprintf(fd,
210 "%ld in %d passes; early out balance=%d\n",
211 (long)a, pass, npts);
212 return a; /* guess was good enough */
213 }
214
215 aa = (sumx/sum-a)*AMP;
216 if (npts >= 2) { /* guess was too low */
217 am = a;
218 aa = xp + max(0.0, aa);
219 if (aa > ap)
220 aa = (a + ap)/2;
221
222 } else if (npts <= -2) { /* guess was two high */
223 ap = a;
224 aa = xm + min(0.0, aa);
225 if (aa < am)
226 aa = (a + am)/2;
227
228 } else {
229 break; /* got it */
230 }
231
232 if (a == aa) {
233 if (trace)
234 fprintf(fd,
235 "%ld in %d passes; force out balance=%d\n",
236 (long)a, pass, npts);
237 return a;
238 }
239 eps = AFAC*abs(aa - a);
240 *eps_ptr = eps;
241 a = aa;
242 }
243
244 if (((x - xlim) % 2) != 0) { /* even number of points? */
245 if (npts == 0) /* yes, return an average */
246 a = (xp+xm)/2;
247 else if (npts > 0)
248 a = (a+xp)/2;
249 else
250 a = (xm+a)/2;
251
252 } else if (npts != 0) { /* odd number of points */
253 if (npts > 0)
254 a = xp;
255 else
256 a = xm;
257 }
258
259 if (trace)
260 fprintf(fd, "%ld in %d passes\n", (long)a, pass);
261 return a;
262 }
263