annotate tools/op.c @ 0:5242703e91d3 tip

Initial checkin for AIM92 aimR8.2 (last updated May 1997).
author tomwalters
date Fri, 20 May 2011 15:19:45 +0100
parents
children
rev   line source
tomwalters@0 1 /*
tomwalters@0 2 op.c Ordered sequence of operations on input stream.
tomwalters@0 3 ------
tomwalters@0 4
tomwalters@0 5 Operations on each input item are done in the order they appear on the
tomwalters@0 6 command line. Operations can be repeated as many times as needed on the
tomwalters@0 7 command line.
tomwalters@0 8
tomwalters@0 9 Input items may be binary numbers (short or float) or ascii lines
tomwalters@0 10 (terminated by CR) according to the type option.
tomwalters@0 11
tomwalters@0 12 The list of operations is as follows. Some operations take valued
tomwalters@0 13 arguments (which are real numbers), others are flags, as indicated:
tomwalters@0 14
tomwalters@0 15 operation result for each input item x
tomwalters@0 16 ----------------- --------------------------------------------
tomwalters@0 17 add=<value> x+<value>
tomwalters@0 18 negate=on -x
tomwalters@0 19 multiply=<value> x*<value>
tomwalters@0 20 divide=<value> x/<value>
tomwalters@0 21 remainder=<value> real remainder of x/<value>
tomwalters@0 22 inverse=on 1/x
tomwalters@0 23 round=on x rounded to nearest integer
tomwalters@0 24 absolute=on |x|
tomwalters@0 25 power=<value> x^<value>
tomwalters@0 26 exponent=<value> <value>^x
tomwalters@0 27 e-exponent=on e^x
tomwalters@0 28 log=<value> log(x) to base <value>
tomwalters@0 29 e-log=on log(x) to base e.
tomwalters@0 30 sin=on sin(x)
tomwalters@0 31 cos=on cos(x)
tomwalters@0 32 tan=on tan(x)
tomwalters@0 33 threshold=<value> if ( x <value> ) 0 else if ( x = <value> ) 1
tomwalters@0 34 diff1=on first difference
tomwalters@0 35 diff2=on second difference
tomwalters@0 36 diff3=on third difference
tomwalters@0 37 cusum=on cumulative sum
tomwalters@0 38 cumean=on cumulative (recursive) mean an
tomwalters@0 39 cumin=on cumulative min
tomwalters@0 40 cumax=on cumulative max
tomwalters@0 41
tomwalters@0 42 The `type' option sets the data type for all input and output items.
tomwalters@0 43 Recognised data types are: short, float, ascii.
tomwalters@0 44
tomwalters@0 45
tomwalters@0 46 Ascii data is read in line-by-line, and comment lines are allowed.
tomwalters@0 47 The following operation apply only to ascii input items, enabling
tomwalters@0 48 comment lines to be stripped or echoed directly to the output.
tomwalters@0 49
tomwalters@0 50 strip=<value> ignore lines beginning with <value> string
tomwalters@0 51 echo=<value> echo lines beginning with <value> string
tomwalters@0 52 number=on number output lines
tomwalters@0 53 If number=on then stripped lines are
tomwalters@0 54 not numbered; they simply dissappear.
tomwalters@0 55
tomwalters@0 56 If there are no operations, (other than comment stripping/echoing or line
tomwalters@0 57 numbering), then ascii lines are echoed.
tomwalters@0 58
tomwalters@0 59
tomwalters@0 60 Examples:
tomwalters@0 61
tomwalters@0 62 1. To specify the following ordered sequence of operations on each binary
tomwalters@0 63 float in the input stream: invert (take reciprocal), scale up by 100,
tomwalters@0 64 square (raise to power 2), invert (take reciprocal), square root (raise to
tomwalters@0 65 power 0.5).
tomwalters@0 66
tomwalters@0 67 op type=float inverse=on multiply=100 power=2 inverse=on power=.5 file
tomwalters@0 68
tomwalters@0 69
tomwalters@0 70 2. To use as an interactive "desk calculator", eg to print out the result
tomwalters@0 71 of log2(17/16):
tomwalters@0 72
tomwalters@0 73 echo 17 | op type=ascii -div16 -log2
tomwalters@0 74
tomwalters@0 75 3. To interactively operate on a list, eg add 1 to each of the list 1,2,3,
tomwalters@0 76 (note how to put newlines into echo'd lists):
tomwalters@0 77
tomwalters@0 78 echo "1\
tomwalters@0 79 2\
tomwalters@0 80 3" | op type=ascii -add1
tomwalters@0 81
tomwalters@0 82 4. To read ascii data, echoing any comment lines beginning with the string
tomwalters@0 83 "**", and for each number x in infile, o/p the number 149-x, do:
tomwalters@0 84
tomwalters@0 85 op type=ascii echo="**" neg=on -add149 infile > outfile
tomwalters@0 86
tomwalters@0 87 */
tomwalters@0 88
tomwalters@0 89 #include <stdio.h>
tomwalters@0 90 #include <math.h>
tomwalters@0 91 #include "options.h"
tomwalters@0 92 #include "strmatch.h"
tomwalters@0 93
tomwalters@0 94 char applic[] = "Ordered sequence of operations on each item in input stream.\n (Input and output in binary shorts)." ;
tomwalters@0 95
tomwalters@0 96 static char *helpstr, *debugstr,
tomwalters@0 97 *addstr ,
tomwalters@0 98 *negstr ,
tomwalters@0 99 *mulstr ,
tomwalters@0 100 *divstr ,
tomwalters@0 101 *remstr ,
tomwalters@0 102 *invstr ,
tomwalters@0 103 *rndstr ,
tomwalters@0 104 *absstr ,
tomwalters@0 105 *powstr ,
tomwalters@0 106 *expstr ,
tomwalters@0 107 *eexpstr,
tomwalters@0 108 *logstr ,
tomwalters@0 109 *elogstr,
tomwalters@0 110 *sinstr ,
tomwalters@0 111 *cosstr ,
tomwalters@0 112 *tanstr ,
tomwalters@0 113 *threshstr,
tomwalters@0 114 *diff1str ,
tomwalters@0 115 *diff2str ,
tomwalters@0 116 *diff3str ,
tomwalters@0 117 *cusumstr ,
tomwalters@0 118 *cumeanstr ,
tomwalters@0 119 *cuminstr ,
tomwalters@0 120 *cumaxstr ,
tomwalters@0 121 *typestr ,
tomwalters@0 122 *stripstr ,
tomwalters@0 123 *echostr ,
tomwalters@0 124 *numstr ;
tomwalters@0 125
tomwalters@0 126 static Options option[] = {
tomwalters@0 127 { "help" , "off" , &helpstr , "help" , DEBUG },
tomwalters@0 128 { "debug" , "off" , &debugstr , "debugging switch" , DEBUG },
tomwalters@0 129 { "add" , "off" , &addstr , "Add x = x+value " , VAL },
tomwalters@0 130 { "negate" , "off" , &negstr , "Negate x = -x " , SETFLAG },
tomwalters@0 131 { "multiply" , "off" , &mulstr , "Multiply x = x*value " , VAL },
tomwalters@0 132 { "divide" , "off" , &divstr , "Divide x = x/value " , VAL },
tomwalters@0 133 { "remainder" , "off" , &remstr , "Remainder x = real remainder of x/value" , VAL },
tomwalters@0 134 { "inverse" , "off" , &invstr , "Inverse x = 1/x " , SETFLAG },
tomwalters@0 135 { "round" , "off" , &rndstr , "Round x = x, to nearest integer" , SETFLAG },
tomwalters@0 136 { "absolute" , "off" , &absstr , "Absolute x = |x| " , SETFLAG },
tomwalters@0 137 { "power" , "off" , &powstr , "Power x = x^value " , VAL },
tomwalters@0 138 { "exponent" , "off" , &expstr , "Exponent x = value^x " , VAL },
tomwalters@0 139 { "e-exponent", "off" , &eexpstr , "e-exponent x = e^x " , SETFLAG },
tomwalters@0 140 { "log" , "off" , &logstr , "Log x = log(x), base `value'." , VAL },
tomwalters@0 141 { "e-log" , "off" , &elogstr , "elog x = log(x), base e. " , SETFLAG },
tomwalters@0 142 { "sin" , "off" , &sinstr , "sin x = sin(x) " , SETFLAG },
tomwalters@0 143 { "cos" , "off" , &cosstr , "cos x = cos(x) " , SETFLAG },
tomwalters@0 144 { "tan" , "off" , &tanstr , "tan x = tan(x) " , SETFLAG },
tomwalters@0 145 { "threshold" , "off" , &threshstr , "Threshold x = 0 (x<value) or 1 (x>=value)", VAL },
tomwalters@0 146 { "diff1" , "off" , &diff1str , "First difference " , SETFLAG },
tomwalters@0 147 { "diff2" , "off" , &diff1str , "Second difference " , SETFLAG },
tomwalters@0 148 { "diff3" , "off" , &diff1str , "Third difference " , SETFLAG },
tomwalters@0 149 { "cusum" , "off" , &cusumstr , "Cumulative sum x = current sum" , SETFLAG },
tomwalters@0 150 { "cumean" , "off" , &cumeanstr , "Cumulative mean x = current (recursive) mean" , SETFLAG },
tomwalters@0 151 { "cumin" , "off" , &cuminstr , "Cumulative min x = current min " , SETFLAG },
tomwalters@0 152 { "cumax" , "off" , &cumaxstr , "Cumulative max x = current max " , SETFLAG },
tomwalters@0 153 { "type" , "short" , &typestr , "Data type (short, float, ascii)" , VAL },
tomwalters@0 154 { "strip" , "off" , &stripstr , "Ignore lines beginning with value string" , VAL },
tomwalters@0 155 { "echo" , "off" , &echostr , "Echo lines beginning with value string" , VAL },
tomwalters@0 156 { "number" , "off" , &numstr , "Number output lines" , SETFLAG },
tomwalters@0 157 ( char * ) 0 } ;
tomwalters@0 158
tomwalters@0 159
tomwalters@0 160 #define SIZE 64 /* Max number of program command-line arguments */
tomwalters@0 161
tomwalters@0 162 #define Add 2 /* Functions */
tomwalters@0 163 #define Negate 3
tomwalters@0 164 #define Multiply 4
tomwalters@0 165 #define Divide 5
tomwalters@0 166 #define Remainder 6
tomwalters@0 167 #define Inverse 7
tomwalters@0 168 #define Round 8
tomwalters@0 169 #define Absolute 9
tomwalters@0 170 #define Power 10
tomwalters@0 171 #define Exp 11
tomwalters@0 172 #define Eexp 12
tomwalters@0 173 #define Log 13
tomwalters@0 174 #define Elog 14
tomwalters@0 175 #define Sin 15
tomwalters@0 176 #define Cos 16
tomwalters@0 177 #define Tan 17
tomwalters@0 178 #define Threshold 18
tomwalters@0 179 #define Diff1 19
tomwalters@0 180 #define Diff2 20
tomwalters@0 181 #define Diff3 21
tomwalters@0 182 #define Sum 22
tomwalters@0 183 #define Mean 23
tomwalters@0 184 #define Min 24
tomwalters@0 185 #define Max 25
tomwalters@0 186 #define Ascii 26
tomwalters@0 187 #define Strip 27
tomwalters@0 188 #define Echo 28
tomwalters@0 189 #define Number 29
tomwalters@0 190
tomwalters@0 191
tomwalters@0 192 /* data types */
tomwalters@0 193
tomwalters@0 194 #define SHORT 0
tomwalters@0 195 #define FLOAT 1
tomwalters@0 196 #define ASCII 2
tomwalters@0 197
tomwalters@0 198 int Type ;
tomwalters@0 199
tomwalters@0 200 int STRIP=0;
tomwalters@0 201 char Stripstr[]="*"; /* Default strip comment string */
tomwalters@0 202
tomwalters@0 203 int ECHO=0;
tomwalters@0 204 char Echostr[]="**"; /* Default echo comment string */
tomwalters@0 205
tomwalters@0 206 int NUMBER=0; /* Flag for numbering ascii lines */
tomwalters@0 207 int n=0;
tomwalters@0 208
tomwalters@0 209 /*
tomwalters@0 210 Rounding for both +ve and -ve numbers. Real numbers are truncated to ints in
tomwalters@0 211 the direction of zero (EG both 0.9 and -0.9 are truncated to 0). Using round,
tomwalters@0 212 0.9 becomes 1 and -0.9 becomes -1. (Note also 0.5 becomes 1, -0.5 becomes -1.
tomwalters@0 213 */
tomwalters@0 214 #define round(x) ((x >= 0) ? ((int)(x+0.5)) : ((int)(x-0.5)))
tomwalters@0 215
tomwalters@0 216 /* Remainder of real division for +ve and -ve numbers */
tomwalters@0 217 #define rem(x,y) ((x >= 0) ? (fabs((x) - (y)*(int)((x)/(y)))) : (fabs((x) - (y)*(int)((x)/(y)-0.5))))
tomwalters@0 218
tomwalters@0 219 /* Threshold and convert to binary (0,1) value */
tomwalters@0 220 #define threshold(x,T) ((x >= T) ? 1 : 0)
tomwalters@0 221
tomwalters@0 222 float diff1();
tomwalters@0 223 float diff2();
tomwalters@0 224 float diff3();
tomwalters@0 225 float mean();
tomwalters@0 226 float min();
tomwalters@0 227 float max();
tomwalters@0 228 float sum();
tomwalters@0 229
tomwalters@0 230 struct opstr { /* array of operations: functions and (optional) arg */
tomwalters@0 231 char func;
tomwalters@0 232 float arg;
tomwalters@0 233 } op[SIZE];
tomwalters@0 234
tomwalters@0 235 int k=0; /* number of operations stored in array */
tomwalters@0 236
tomwalters@0 237 main(argc, argv)
tomwalters@0 238 int argc ;
tomwalters@0 239 char *argv[] ;
tomwalters@0 240 {
tomwalters@0 241 FILE *fp, *fopen();
tomwalters@0 242 float xf, ops();
tomwalters@0 243 short xs ;
tomwalters@0 244 char s[256], stripcomment[32], echocomment[32];
tomwalters@0 245 char *prog, *val ;
tomwalters@0 246 int i, index, span ;
tomwalters@0 247
tomwalters@0 248 strcpy(stripcomment, Stripstr);
tomwalters@0 249 strcpy(echocomment, Echostr );
tomwalters@0 250
tomwalters@0 251 for (i=0 ; option[i].name != (char *)0 ; i++)
tomwalters@0 252 *(option[i].val) = option[i].dflt ;
tomwalters@0 253
tomwalters@0 254 prog = argv[0] ;
tomwalters@0 255 while ( --argc > 0 && ++argv ) {
tomwalters@0 256 if ( ( i = whichopt( option, *argv, &index, &span ) ) == UNKNOWN )
tomwalters@0 257 break ;
tomwalters@0 258 if ( ( val = checksyntax( *argv, span, option[index].type ) ) == (char *)0 )
tomwalters@0 259 break ;
tomwalters@0 260 if ( i == AMBIGUOUS ) {
tomwalters@0 261 fprintf(stderr,"%s: ambiguous option [%s]\n", prog, *argv ) ;
tomwalters@0 262 exit ( 1 ) ;
tomwalters@0 263 }
tomwalters@0 264 switch ( index ) {
tomwalters@0 265 case 0 : operate( option, index, val ) ; break;
tomwalters@0 266 case 1 : operate( option, index, val ) ; break;
tomwalters@0 267 case 2 : op[k].func = Add; op[k].arg = atof(val); k++; break;
tomwalters@0 268 case 3 : op[k].func = Negate; k++; break;
tomwalters@0 269 case 4 : op[k].func = Multiply; op[k].arg = atof(val); k++; break;
tomwalters@0 270 case 5 : op[k].func = Divide; op[k].arg = atof(val); k++; break;
tomwalters@0 271 case 6 : op[k].func = Remainder; op[k].arg = atof(val); k++; break;
tomwalters@0 272 case 7 : op[k].func = Inverse; k++; break;
tomwalters@0 273 case 8 : op[k].func = Round; k++; break;
tomwalters@0 274 case 9 : op[k].func = Absolute; k++; break;
tomwalters@0 275 case 10 : op[k].func = Power; op[k].arg = atof(val); k++; break;
tomwalters@0 276 case 11 : op[k].func = Exp; op[k].arg = atof(val); k++; break;
tomwalters@0 277 case 12 : op[k].func = Eexp; k++; break;
tomwalters@0 278 case 13 : op[k].func = Log; op[k].arg = atof(val); k++; break;
tomwalters@0 279 case 14 : op[k].func = Elog; k++; break;
tomwalters@0 280 case 15 : op[k].func = Sin; k++; break;
tomwalters@0 281 case 16 : op[k].func = Cos; k++; break;
tomwalters@0 282 case 17 : op[k].func = Tan; k++; break;
tomwalters@0 283 case 18 : op[k].func = Threshold; op[k].arg = atof(val); k++; break;
tomwalters@0 284 case 19 : op[k].func = Diff1; k++; break;
tomwalters@0 285 case 20 : op[k].func = Diff2; k++; break;
tomwalters@0 286 case 21 : op[k].func = Diff3; k++; break;
tomwalters@0 287 case 22 : op[k].func = Sum; k++; break;
tomwalters@0 288 case 23 : op[k].func = Mean; k++; break;
tomwalters@0 289 case 24 : op[k].func = Min; k++; break;
tomwalters@0 290 case 25 : op[k].func = Max; k++; break;
tomwalters@0 291 case 26 : Type = checktype( val ) ; break;
tomwalters@0 292 case 27 : STRIP++; strcpy(stripcomment, val); break;
tomwalters@0 293 case 28 : ECHO++; strcpy(echocomment, val); break;
tomwalters@0 294 case 29 : NUMBER++; break;
tomwalters@0 295 }
tomwalters@0 296
tomwalters@0 297 }
tomwalters@0 298 if ( !isoff( helpstr ) )
tomwalters@0 299 helpopts( helpstr, prog, applic, option ) ;
tomwalters@0 300
tomwalters@0 301 if (argc <= 0) fp = stdin;
tomwalters@0 302 else if ((fp = fopen(*argv, "r")) == NULL) {
tomwalters@0 303 fprintf(stderr,"can't open %s\n", *argv);
tomwalters@0 304 exit(1);
tomwalters@0 305 }
tomwalters@0 306
tomwalters@0 307 switch ( Type ) {
tomwalters@0 308
tomwalters@0 309 case SHORT :
tomwalters@0 310 while ( fread( &xs, sizeof(short), 1, fp ) ) {
tomwalters@0 311 xs = (short)ops( (float)xs ) ;
tomwalters@0 312 fwrite( &xs, sizeof(short), 1, stdout ) ;
tomwalters@0 313 }
tomwalters@0 314 break ;
tomwalters@0 315
tomwalters@0 316 case FLOAT :
tomwalters@0 317 while ( fread( &xf, sizeof(float), 1, fp ) ) {
tomwalters@0 318 xf = ops( xf ) ;
tomwalters@0 319 fwrite( &xf, sizeof(float), 1, stdout ) ;
tomwalters@0 320 }
tomwalters@0 321 break ;
tomwalters@0 322
tomwalters@0 323 case ASCII :
tomwalters@0 324 while (fgets(s, 256, fp)) {
tomwalters@0 325 if (STRIP && match(s,stripcomment)==0)
tomwalters@0 326 ; /* strip comment (ie do nothing) */
tomwalters@0 327 else {
tomwalters@0 328 if (NUMBER) printf("%d: ", ++n);
tomwalters@0 329 if (ECHO && match(s,echocomment)==0)
tomwalters@0 330 fputs(s,stdout); /* echo comment */
tomwalters@0 331 else {
tomwalters@0 332 if ( k>0 ) {
tomwalters@0 333 xf = ops( atof(s) ) ;
tomwalters@0 334 printf("%.3f\n", xf ) ;
tomwalters@0 335 }
tomwalters@0 336 else fputs(s,stdout); /* echo line (case of no ops) */
tomwalters@0 337 }
tomwalters@0 338 }
tomwalters@0 339 }
tomwalters@0 340 break ;
tomwalters@0 341 }
tomwalters@0 342
tomwalters@0 343 fclose(fp);
tomwalters@0 344 }
tomwalters@0 345
tomwalters@0 346
tomwalters@0 347 checktype( s )
tomwalters@0 348 char *s ;
tomwalters@0 349 {
tomwalters@0 350 if ( iststr( s, "short" ) ) return SHORT ;
tomwalters@0 351 if ( iststr( s, "float" ) ) return FLOAT ;
tomwalters@0 352 if ( iststr( s, "ascii" ) ) return ASCII ;
tomwalters@0 353 fprintf( stderr,"unknown datatype [%s]\n", s ) ;
tomwalters@0 354 exit( 1 ) ;
tomwalters@0 355 }
tomwalters@0 356
tomwalters@0 357
tomwalters@0 358 float ops(x)
tomwalters@0 359 float x;
tomwalters@0 360 {
tomwalters@0 361 int i;
tomwalters@0 362
tomwalters@0 363 for (i=0 ; i<k ; i++)
tomwalters@0 364 switch (op[i].func) {
tomwalters@0 365 case Add: x = x+op[i].arg; break;
tomwalters@0 366 case Negate: x = (-x); break;
tomwalters@0 367 case Multiply: x = x*op[i].arg; break;
tomwalters@0 368 case Divide: x = x/op[i].arg; break;
tomwalters@0 369 case Remainder: x = rem(x,op[i].arg); break;
tomwalters@0 370 case Inverse: x = 1.0/x; break;
tomwalters@0 371 case Round: x = round(x); break;
tomwalters@0 372 case Absolute: x = fabs(x); break;
tomwalters@0 373 case Power: x = pow(x,op[i].arg); break;
tomwalters@0 374 case Exp: x = pow(op[i].arg,x); break;
tomwalters@0 375 case Eexp: x = exp(x); break;
tomwalters@0 376 case Log: x = log(x)/log(op[i].arg); break;
tomwalters@0 377 case Elog: x = log(x); break;
tomwalters@0 378 case Sin: x = sin(x); break;
tomwalters@0 379 case Cos: x = cos(x); break;
tomwalters@0 380 case Tan: x = tan(x); break;
tomwalters@0 381 case Threshold: x = threshold(x,op[i].arg); break;
tomwalters@0 382 case Diff1: x = diff1(x); break;
tomwalters@0 383 case Diff2: x = diff2(x); break;
tomwalters@0 384 case Diff3: x = diff3(x); break;
tomwalters@0 385 case Sum: x = sum(x); break;
tomwalters@0 386 case Mean: x = mean(x); break;
tomwalters@0 387 case Min: x = min(x); break;
tomwalters@0 388 case Max: x = max(x); break;
tomwalters@0 389 }
tomwalters@0 390 return x;
tomwalters@0 391 }
tomwalters@0 392
tomwalters@0 393
tomwalters@0 394 match(s1,s2)
tomwalters@0 395 char *s1, *s2;
tomwalters@0 396 {
tomwalters@0 397 int n1, n2;
tomwalters@0 398
tomwalters@0 399 if ((n1=strlen(s1)) < (n2=strlen(s2)))
tomwalters@0 400 return strncmp(s1,s2,n1);
tomwalters@0 401 else
tomwalters@0 402 return strncmp(s1,s2,n2);
tomwalters@0 403 }
tomwalters@0 404
tomwalters@0 405
tomwalters@0 406 /****************************************************************************
tomwalters@0 407 The past values used to compute differences are, (in the order they would
tomwalters@0 408 appear in an array, where x is the current value):
tomwalters@0 409 1st diffs ... x11 x
tomwalters@0 410 2nd diffs ... x22 x21 x
tomwalters@0 411 3rd diffs ... x33 x32 x31 x
tomwalters@0 412 *****************************************************************************/
tomwalters@0 413 float diff1(x)
tomwalters@0 414 float x;
tomwalters@0 415 {
tomwalters@0 416 static float x11=0;
tomwalters@0 417 float d;
tomwalters@0 418
tomwalters@0 419 d=fabs(x11-x);
tomwalters@0 420 x11=x;
tomwalters@0 421 return d;
tomwalters@0 422 }
tomwalters@0 423
tomwalters@0 424 float diff2(x)
tomwalters@0 425 float x;
tomwalters@0 426 {
tomwalters@0 427 static float x21=0, x22=0;
tomwalters@0 428 float d;
tomwalters@0 429
tomwalters@0 430 d=fabs( fabs(x22-x21)-fabs(x21-x) );
tomwalters@0 431 x22=x21;
tomwalters@0 432 x21=x;
tomwalters@0 433 return d;
tomwalters@0 434 }
tomwalters@0 435
tomwalters@0 436 float diff3(x)
tomwalters@0 437 float x;
tomwalters@0 438 {
tomwalters@0 439 static float x31=0, x32=0, x33=0;
tomwalters@0 440 float d;
tomwalters@0 441
tomwalters@0 442 d=fabs( fabs( fabs(x33-x32)-fabs(x32-x31) ) - fabs( fabs(x32-x31)-fabs(x31-x) ) );
tomwalters@0 443 x33=x32;
tomwalters@0 444 x32=x31;
tomwalters@0 445 x31=x;
tomwalters@0 446 return d;
tomwalters@0 447 }
tomwalters@0 448
tomwalters@0 449 /****************************************************************************
tomwalters@0 450 Recurisve estimation of the mean
tomwalters@0 451 (Ref: Young, eqn 2.14 (p14), eqn VI(2) (p63), eqn 5.17 (p65))
tomwalters@0 452 ****************************************************************************/
tomwalters@0 453 float mean(x)
tomwalters@0 454 float x;
tomwalters@0 455 {
tomwalters@0 456 static float a = 0; /* initial mean value */
tomwalters@0 457 static float p = 10000; /* initial gain factor [Young p65] */
tomwalters@0 458
tomwalters@0 459 p = p/(1+p); /* update p (the gain factor) */
tomwalters@0 460 a = a - p*(a-x); /* update a (the mean estimate) */
tomwalters@0 461 return a;
tomwalters@0 462 }
tomwalters@0 463
tomwalters@0 464 /****************************************************************************
tomwalters@0 465 Min, Max, and Sum
tomwalters@0 466 ****************************************************************************/
tomwalters@0 467 float min(x)
tomwalters@0 468 float x;
tomwalters@0 469 {
tomwalters@0 470 static float MIN=1.0e30;
tomwalters@0 471
tomwalters@0 472 if (x < MIN) MIN = x;
tomwalters@0 473 return MIN;
tomwalters@0 474 }
tomwalters@0 475
tomwalters@0 476 float max(x)
tomwalters@0 477 float x;
tomwalters@0 478 {
tomwalters@0 479 static float MAX=0;
tomwalters@0 480
tomwalters@0 481 if (x > MAX) MAX = x;
tomwalters@0 482 return MAX;
tomwalters@0 483 }
tomwalters@0 484
tomwalters@0 485 float sum(x)
tomwalters@0 486 float x;
tomwalters@0 487 {
tomwalters@0 488 static float SUM=0;
tomwalters@0 489
tomwalters@0 490 SUM += x;
tomwalters@0 491 return SUM;
tomwalters@0 492 }