/* GRASP: Copyright 1997,1998 Bruce Allen */ static char *rcsid="\$Id: utility.c,v 1.23 1999/07/01 20:51:41 ballen Exp \$\ \n\$Name: \$"; #include "grasp.h" #ifndef _WIN32 #include #endif /* when called with pointers to arrays c,a,b, this routine takes sets c = a x b where a,b,c are complex numbers */ void product(float *c,float *a, float *b,int ncomplex) { float ar,ai,br,bi; while (ncomplex-->0) { ar=*a++; ai=*a++; br=*b++; bi=*b++; *c++=ar*br-ai*bi; *c++=ar*bi+ai*br; } return; } /* when called with pointers to arrays c,a,b, this routine sets c = a / b where a,b,c are complex numbers */ void ratio(float *c,float *a, float *b,int ncomplex) { float ar,ai,br,bi,mod; while (ncomplex-->0) { ar=*a++; ai=*a++; br=*b++; bi=*b++; mod=1.0/(br*br+bi*bi); *c++=mod*(ar*br+ai*bi); *c++=mod*(ai*br-ar*bi); } return; } /* when called with pointers to arrays c,a,b, this routine sets c = a x b* where a,b,c are complex numbers, and * is conjugate. */ void productc(float *c,float *a, float *b,int ncomplex) { float ar,ai,br,bi; while (ncomplex-->0) { ar=*a++; ai=*a++; br=*b++; bi=*b++; *c++=ar*br+ai*bi; *c++=ai*br-ar*bi; } return; } void avg_spec(float *data, float *average, int npoint, int *reset, float srate, float decaytime, int windowtype, int overlap) { void realft(float*,unsigned long,int); static float *saved=NULL,*working=NULL,*window=NULL; static double wss,decay,lastnorm; float spec,*dataptr; double win,delta,norm; int i,ipn,im,re; /* Keep the data */ dataptr=data; /* if we are resetting, then.. */ if (*reset==0) { /* allocate space for saved data */ saved=(float *)realloc(saved,sizeof(float)*npoint/2); if (saved==NULL) { GR_start_error("avg_spec()",rcsid,__FILE__,__LINE__); GR_report_error("failed to allocate %d floats\n",npoint); GR_end_error(); abort(); } /* allocate working space */ working=(float *)realloc(working,sizeof(float)*npoint); if (working==NULL) { GR_start_error("avg_spec()",rcsid,__FILE__,__LINE__); GR_report_error("failed to allocate %d floats\n",npoint); GR_end_error(); abort(); } /* allocate space for window function */ window=(float *)realloc(window,sizeof(float)*npoint); if (window==NULL) { GR_start_error("avg_spec()",rcsid,__FILE__,__LINE__); GR_report_error("failed to allocate %d floats\n",npoint); GR_end_error(); abort(); } /* set average to zero and put data in saved space */ for (i=0;i=0){ /* window the data data */ for (i=0;i=nbins) (*over)++; else /* and if it's in range, then bin it! */ bins[offset]+=1.0; } return; } /* Here input points to an array of shorts, of length ninput. This routines takes these values, input[0..ninput-1] and bins them into bins. Each element input[r] increments the value of bins[input[r]+offset] by one. */ void binshort(short *input,int ninput,double *bins,int offset){ int i; /* step through the input array */ for (i=0;i/dev/null 2>&1 &"); #endif /* and return */ return; } void graph_double(double *array,int n,int spacing) { FILE *fp; int i; /* open a file for writing */ fp=fopen("temp.graph","w"); if (fp==NULL) { GR_start_error("graph_double()",rcsid,__FILE__,__LINE__); GR_report_error("Can't write file \"temp.graph\" in current directory -- sorry!\n"); GR_end_error(); return; } /* print data into the file */ for (i=0;i/dev/null 2>&1 &"); #endif /* and return */ return; } void graph_short(short *array,int n) { FILE *fp; int i; /* open a file for writing */ fp=fopen("temp.graph","w"); if (fp==NULL) { GR_start_error("graph_short()",rcsid,__FILE__,__LINE__); GR_report_error("Can't write file \"temp.graph\" in current directory -- sorry!\n"); GR_end_error(); return; } /* print data into the file */ for (i=0;i/dev/null 2>&1 &"); #endif /* and return */ return; } void sgraph(short *array,int n,char* ext,int filenumber) { FILE *fp; int i; char name[256]; /* open a file for writing */ sprintf(name,"%s.%03d",ext,filenumber); fp=fopen(name,"w"); if (fp==NULL) { GR_start_error("sgraph()",rcsid,__FILE__,__LINE__); GR_report_error("Can't write file %s in current directory -- sorry!\n",name); GR_end_error(); return; } /* print data into the file */ for (i=0;istorage+3*lastbins) */ histogram=storage+2*lastbins; /* clear out the histogram */ for (i=0;iexceed+2.0*sqrt(exceed) || sum01+sum13+sum35!=n) retval=0; /* print info if desired */ if (print) printf(" Distribution: s= %d, N>3s= %d (expect %d), N>5s= %d (expect 0)\n", sigma,(int)sum35,(int)exceed,(int)(n-sum01-sum35-sum13)); return retval; } /* This routine opens files, appending the path determined by the environment variable with the character string shortpath. It tries to print intelligent error messages if the environment variable is not defined, or if the file pointed to does not exist. */ FILE* grasp_open(const char *environment_variable,const char *shortpath,const char *mode) { char *pathhead; char longpath[256]; FILE *fp; /* get the environment variable giving the data path */ pathhead=getenv(environment_variable); /* if this is null, then print out a warning to the user */ if (!pathhead) { GR_start_error("grasp_open()",rcsid,__FILE__,__LINE__); GR_report_error("the environment variable %s must be set.\n",environment_variable); GR_report_error("It should point to a directory containing 40-meter data, or parameters.\n"); GR_report_error("It can be set with a command like:\nsetenv %s /this/is/the/path\n",environment_variable); GR_end_error(); exit(1); } /* now copy this (head of the data path) into the return */ strcpy(longpath,pathhead); /* and append the remainder of the path */ #ifdef __MACOS__ strcat(longpath,":"); #elif defined _WIN32 strcat(longpath,"\\"); #else strcat(longpath,"/"); #endif strcat(longpath,shortpath); /* now open the file */ fp=fopen(longpath,mode); if (fp==NULL) { GR_start_error("grasp_open()",rcsid,__FILE__,__LINE__); GR_report_error("Unable to open file %s in mode \"%s\"\n",longpath,mode); GR_report_error("The current value of environment variable %s is %s\n",environment_variable,pathhead); GR_end_error(); exit(1); } return fp; }