00001 #include "ctypes_def.h"
00002 #if TYPE == FLOAT
00003 #define REJECT sreject
00004 #elif TYPE == DOUBLE
00005 #define REJECT dreject
00006 #elif TYPE == COMPLEXFLOAT
00007 #define REJECT creject
00008 #elif TYPE == COMPLEXDOUBLE
00009 #define REJECT zreject
00010 #endif
00011
00012
00013 void REJECT(int n, CTYPE *data, int *isgood, RTYPE factor)
00014 {
00015 #ifndef NDEBUG
00016 int verbose=1;
00017 #else
00018 int verbose=0;
00019 #endif
00020 int i,ngood;
00021 CTYPE sum;
00022 RTYPE std,sum2;
00023
00024 ngood = 0;
00025 sum=CZERO;
00026 sum2=RZERO;
00027 for (i=0; i<n; i++)
00028 {
00029 if (isgood[i])
00030 {
00031 sum += data[i];
00032 sum2 += CONJ(data[i])*data[i];
00033 ++ngood;
00034 }
00035 }
00036 std = sqrt(fabs((sum2-CONJ(sum)*sum/ngood)/(ngood-1)));
00037 if (verbose)
00038 printf("ngood = %d, sum2 = %f, sum^2 = %f, std = %f\n",
00039 ngood,sum2,SQR(sum),std);
00040
00041
00042 for (i=0; i<n; i++)
00043 {
00044 if (isgood[i] && fabs(data[i])>factor*std)
00045 {
00046 data[i] = CZERO;
00047 isgood[i] = 0;
00048 }
00049 }
00050 }
00051
00052 #undef REJECT
00053 #include "ctypes_undef.h"
