home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1994 January / usenetsourcesnewsgroupsinfomagicjanuary1994.iso / sources / misc / volume21 / cloops / part01 / statw.c < prev    next >
Encoding:
C/C++ Source or Header  |  1991-07-25  |  2.1 KB  |  87 lines
  1. /*
  2.  * This file is part of the Livermore Loops transliteration into C.
  3.  * Copyright (C) 1991 by Martin Fouts
  4.  *
  5.  * This program is free software; you can redistribute it and/or modify
  6.  * it under the terms of the GNU General Public License as published by
  7.  * the Free Software Foundation; either version 1, or (at your option)
  8.  * any later version.
  9.  *
  10.  * This program is distributed in the hope that it will be useful,
  11.  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  12.  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13.  * GNU General Public License for more details.
  14.  *
  15.  * You should have received a copy of the GNU General Public License
  16.  * along with this program; if not, write to the Free Software
  17.  * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  18.  */
  19.  
  20. #include "types.h"
  21.  
  22. Void sordid();
  23. Float sqrt();
  24.  
  25. Void statw(result,ox,ix,x,w,n)
  26. Float result[], ox[], x[], w[];
  27. Int ix[], n;
  28. {
  29.   Int k;
  30.   Float a, r, s, t;
  31.   Float d, u, v, h;
  32.   for (k = 0; k < 9; k++)
  33.     result[k] = 0.0;
  34.   if (n <= 0) return;
  35.   a = 0.0; s = 0.0; t = 0.0;
  36.   for (k = 0; k < n; k++) {
  37.     s += w[k] * x[k];
  38.     t += w[k];
  39.   }
  40.   if (t != 0.0) a = s / t;
  41.   result[0] = a;
  42.   d = 0.0; u = 0.0;
  43.   for (k = 0; k < n; k++)
  44.     d += w[k] * (x[k]-a) * (x[k] - a);
  45.   if (t != 0.0) d = d/t;
  46.   if (d >= 0.0) u = sqrt(d);
  47.   result[1] = u;
  48.   u = x[0]; v = x[0];
  49.   for (k = 1; k < n; k++) {
  50.     if (u > x[k]) u = x[k];
  51.     if (v < x[k]) v = x[k];
  52.   }
  53.   result[2] = u;
  54.   result[3] = v;
  55.   h = 0.0;
  56.   for (k = 0; k < n; k++)
  57.     if (x[k] != 0.0)
  58.       h += w[k]/x[k];
  59.   if (h != 0.0) h = t / h;
  60.   result[4] = h;
  61.   result[5] = t;
  62.   sordid( ix, ox, x, n, (Int)1);
  63.   r = 0.0;
  64.   for (k = 0; k < n; k++) {
  65.     r += w[ix[k]];
  66.     if (r > 0.5 * t) break;
  67.   }
  68.   result[6] = ox[k];
  69.   result[7] = (Float)k;
  70.   for (k = 0; k < n; k++) {
  71.     Float temp;
  72.     temp = x[k] - result[6];
  73.     if (temp < 0.0) temp = -temp;
  74.     ox[k] = temp;
  75.   }
  76.   sordid(ix, ox, ox, n, (Int)1);
  77.   r = 0.0;
  78.   for (k = 0; k < n; k++) {
  79.     r += w[ix[k]];
  80.     if (r > 0.7 * t) break;
  81.   }
  82.   result[8] = ox[k];
  83.   sordid(ix, ox, x, n, (Int)2);
  84.   return;
  85. }
  86.  
  87.