usr.sbin/timed/timed/networkdelta.c

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