GRASS Programmer's Manual
6.4.2(2012)
|
00001 #include <stdio.h> 00002 #include <stdlib.h> 00003 #include <math.h> 00004 #include "local_proto.h" 00005 00006 00007 /*- 00008 * driver program for AS 181: Royston's extension of the Shapiro-Wilk 00009 * W statistic to n=2000 00010 * needs as181.c as177.c as241.c dcmp.c as66.c 00011 */ 00012 00013 double *royston(double *x, int n) 00014 { 00015 static double y[2]; 00016 double *a, eps, w, pw, mean = 0, ssq = 0, *xcopy; 00017 int i, ifault, n2; 00018 00019 n2 = (int)floor((double)n / 2); 00020 00021 #ifndef lint 00022 if ((a = (double *)malloc(n2 * sizeof(double))) == NULL) { 00023 fprintf(stderr, "Memory error in royston\n"); 00024 exit(EXIT_FAILURE); 00025 } 00026 if ((xcopy = (double *)malloc(n * sizeof(double))) == NULL) { 00027 fprintf(stderr, "Memory error in royston\n"); 00028 exit(EXIT_FAILURE); 00029 } 00030 #endif /* lint */ 00031 00032 for (i = 0; i < n; ++i) { 00033 xcopy[i] = x[i]; 00034 mean += x[i]; 00035 } 00036 mean /= n; 00037 00038 qsort(xcopy, n, sizeof(double), dcmp); 00039 00040 for (i = 0; i < n; ++i) 00041 ssq += (mean - x[i]) * (mean - x[i]); 00042 00043 wcoef(a, n, n2, &eps, &ifault); 00044 00045 if (ifault == 0) 00046 wext(xcopy, n, ssq, a, n2, eps, &w, &pw, &ifault); 00047 else { 00048 fprintf(stderr, "Error in wcoef()\n"); 00049 return (double *)NULL; 00050 } 00051 00052 if (ifault == 0) { 00053 y[0] = w; 00054 y[1] = pw; 00055 } 00056 else { 00057 fprintf(stderr, "Error in wcoef()\n"); 00058 return (double *)NULL; 00059 } 00060 00061 free(a); 00062 free(xcopy); 00063 00064 return y; 00065 }