SDDS ToolKit Programs and Libraries for C and Python
All Classes Files Functions Variables Macros Pages
sddsenvelope.c
Go to the documentation of this file.
1/**
2 * @file sddsenvelope.c
3 * @brief Combine data from SDDS pages to create a new file with computed statistics.
4 *
5 * @details
6 * This program processes SDDS (Self Describing Data Set) files, performing statistical computations
7 * such as maximum, minimum, mean, and others across pages of data. The resulting statistics are stored
8 * in an output SDDS file for further analysis.
9 *
10 * @section Usage
11 * ```
12 * sddsenvelope [<inputfile>] [<outputfile>]
13 * [-pipe=[input][,output]]
14 * [-nowarnings]
15 * [-maximum=<column-names>]
16 * [-minimum=<column-names>]
17 * [-cmaximum=<indep-column>,<column-names>]
18 * [-cminimum=<indep-column>,<column-names>]
19 * [-pmaximum=<indep-parameter>,<column-names>]
20 * [-pminimum=<indep-parameter>,<column-names>]
21 * [-largest=<column-names>]
22 * [-signedLargest=<column-names>]
23 * [-mean=<column-names>]
24 * [-sum=<power>,<column-names>]
25 * [-median=<column-names>]
26 * [-decilerange=<column-names>]
27 * [-percentile=<percentage>,<column-names>]
28 * [-standarddeviation=<column-names>]
29 * [-rms=<column-names>]
30 * [-sigma=<column-names>]
31 * [-slope=<indep-parameter>,<column-names>]
32 * [-intercept=<indep-parameter>,<column-names>]
33 * [-wmean=<weightColumn>,<columnNames>]
34 * [-wstandarddeviation=<weightColumn>,<columnNames>]
35 * [-wrms=<weightColumn>,<columnNames>]
36 * [-wsigma=<weightColumn>,<columnNames>]
37 * [-majorOrder=row|column]
38 * ```
39 *
40 * @section Options
41 * | Option | Description |
42 * |------------------------------------|--------------------------------------------------------------------------|
43 * | `-pipe` | Use pipe for input/output. |
44 * | `-nowarnings` | Suppress warnings. |
45 * | `-maximum` | Compute maximum values for specified columns. |
46 * | `-minimum` | Compute minimum values for specified columns. |
47 * | `-cmaximum` | Compute conditional maximum based on an independent column. |
48 * | `-cminimum` | Compute conditional minimum based on an independent column. |
49 * | `-pmaximum` | Parameter-based maximum. |
50 * | `-pminimum` | Parameter-based minimum. |
51 * | `-largest` | Compute largest absolute values. |
52 * | `-signedLargest` | Compute largest signed values. |
53 * | `-mean` | Compute mean values. |
54 * | `-sum` | Compute sum raised to a specified power. |
55 * | `-median` | Compute median values. |
56 * | `-decilerange` | Compute decile range. |
57 * | `-percentile` | Compute specified percentile. |
58 * | `-standarddeviation` | Compute standard deviation. |
59 * | `-rms` | Compute root mean square (RMS) values. |
60 * | `-sigma` | Compute sigma values. |
61 * | `-slope` | Compute slope for linear fit. |
62 * | `-intercept` | Compute intercept for linear fit. |
63 * | `-wmean` | Compute weighted mean. |
64 * | `-wstandarddeviation` | Compute weighted standard deviation. |
65 * | `-wrms` | Compute weighted RMS. |
66 * | `-wsigma` | Compute weighted sigma. |
67 * | `-majorOrder` | Set data ordering for output file. |
68 *
69 * @subsection Incompatibilities
70 * - Only one of the following options may be specified:
71 * - `-mean`
72 * - `-median`
73 * - `-percentile`
74 * - `-standarddeviation`
75 *
76 * @subsection SR Specific Requirements
77 * - `-percentile` requires a percentage value between 0 and 100.
78 *
79 * @authors
80 * M. Borland, C. Saunders, R. Soliday, H. Shang
81 *
82 * @license
83 * This file is distributed under the terms of the Software License Agreement
84 * found in the file LICENSE included with this distribution.
85 */
86
87#include "mdb.h"
88#include "scan.h"
89#include "SDDS.h"
90#include <ctype.h>
91
92/* Enumeration for option types */
93enum option_type {
94 SET_COPY,
95 SET_MAXIMA,
96 SET_MINIMA,
97 SET_MEANS,
98 SET_SDS,
99 SET_RMSS,
100 SET_SUMS,
101 SET_SLOPE,
102 SET_INTERCEPT,
103 SET_PIPE,
104 SET_SIGMAS,
105 SET_MEDIAN,
106 SET_DRANGE,
107 SET_WMEANS,
108 SET_WSDS,
109 SET_WRMSS,
110 SET_WSIGMAS,
111 SET_NOWARNINGS,
112 SET_LARGEST,
113 SET_PERCENTILE,
114 SET_SIGNEDLARGEST,
115 SET_PMAXIMA,
116 SET_PMINIMA,
117 SET_MAJOR_ORDER,
118 SET_EXMM_MEAN,
119 SET_CMAXIMA,
120 SET_CMINIMA,
121 N_OPTIONS
122};
123
124char *option[N_OPTIONS] = {
125 "copy",
126 "maximum",
127 "minimum",
128 "mean",
129 "standarddeviations",
130 "rms",
131 "sum",
132 "slope",
133 "intercept",
134 "pipe",
135 "sigmas",
136 "median",
137 "decilerange",
138 "wmean",
139 "wstandarddeviations",
140 "wrms",
141 "wsigma",
142 "nowarnings",
143 "largest",
144 "percentile",
145 "signedlargest",
146 "pmaximum",
147 "pminimum",
148 "majorOrder",
149 "exmmMean",
150 "cmaximum",
151 "cminimum",
152};
153
154char *optionSuffix[N_OPTIONS] = {
155 "",
156 "Max",
157 "Min",
158 "Mean",
159 "StDev",
160 "Rms",
161 "Sum",
162 "Slope",
163 "Intercept",
164 "",
165 "Sigma",
166 "Median",
167 "DRange",
168 "WMean",
169 "WStDev",
170 "WRms",
171 "WSigma",
172 "",
173 "Largest",
174 "Percentile",
175 "SignedLargest",
176 "PMaximum",
177 "PMinimum",
178 "",
179 "ExmmMean",
180 "CMaximum",
181 "CMinimum",
182};
183
184/* this structure stores a command-line request for statistics computation */
185/* columnName may contain wildcards */
186typedef struct
187{
188 char *columnName;
189 char *weightColumnName;
190 long optionCode, sumPower;
191 double percentile;
192 char *percentileString;
193 char *functionOf;
195
196/* this structure stores data necessary for accessing/creating SDDS columns and
197 * for computing a statistic
198 */
199typedef struct
200{
201 char *sourceColumn, *weightColumn, *resultColumn, *functionOf;
202 long optionCode, resultIndex, sumPower;
203 double percentile;
204 char *percentileString;
205 /* these store intermediate values during processing */
206 void *copy;
207 double *value1, *value2, *value3, *value4;
208 double **array;
209 double *sumWeight;
211
212long addStatRequests(STAT_REQUEST **statRequest, long requests, char **item, long items, long code, double percentile, long power, char *functionOf, long weighted, char *percentileString);
213/*weighted=0, no weighted column; else, weighted statistic, the weight factor is given by
214 weightedColumn*/
215STAT_DEFINITION *compileStatDefinitions(SDDS_DATASET *inTable, long *stats, STAT_REQUEST *request, long requests);
216long setupOutputFile(SDDS_DATASET *outTable, char *output, SDDS_DATASET *inTable, STAT_DEFINITION *stat, long stats, int64_t rows, short columnMajorOrder);
217int compute_mean_exclude_min_max(double *value, double *data, long n);
218
219static char *USAGE = "sddsenvelope [<input>] [<output>] [options]\n"
220 " [-pipe=[input][,output]]\n"
221 " [-nowarnings]\n"
222 " [-maximum=<column-names>]\n"
223 " [-minimum=<column-names>]\n"
224 " [-cmaximum=<indep-column>,<column-names>]\n"
225 " [-cminimum=<indep-column>,<column-names>]\n"
226 " [-pmaximum=<indep-parameter>,<column-names>]\n"
227 " [-pminimum=<indep-parameter>,<column-names>]\n"
228 " [-largest=<column-names>]\n"
229 " [-signedLargest=<column-names>]\n"
230 " [-mean=<column-names>]\n"
231 " [-sum=<power>,<column-names>]\n"
232 " [-median=<column-names>]\n"
233 " [-decilerange=<column-names>]\n"
234 " [-percentile=<percentage>,<column-names>]\n"
235 " [-standarddeviation=<column-names>]\n"
236 " [-rms=<column-names>]\n"
237 " [-sigma=<column-names>]\n"
238 " [-slope=<indep-parameter>,<column-names>]\n"
239 " [-intercept=<indep-parameter>,<column-names>]\n"
240 " [-wmean=<weightColumn>,<columnNames>]\n"
241 " [-wstandarddeviation=<weightColumn>,<columnNames>]\n"
242 " [-wrms=<weightColumn>,<columnNames>]\n"
243 " [-wsigma=<weightColumn>,<columnNames>]\n"
244 " [-majorOrder=row|column]\n"
245 "Options:\n"
246 " -copy=<column-names> Copy specified columns.\n"
247 " -pipe=[input][,output] Use pipe for input/output.\n"
248 " -nowarnings Suppress warnings.\n"
249 " -maximum=<column-names> Compute maximum values.\n"
250 " -minimum=<column-names> Compute minimum values.\n"
251 " -cmaximum=<indep-column>,<column-names> Conditional maximum based on an independent column.\n"
252 " -cminimum=<indep-column>,<column-names> Conditional minimum based on an independent column.\n"
253 " -pmaximum=<indep-parameter>,<column-names> Parameter-based maximum.\n"
254 " -pminimum=<indep-parameter>,<column-names> Parameter-based minimum.\n"
255 " -largest=<column-names> Compute the largest absolute values.\n"
256 " -signedLargest=<column-names> Compute the largest signed values.\n"
257 " -mean=<column-names> Compute mean values.\n"
258 " -sum=<power>,<column-names> Compute sum with power.\n"
259 " -median=<column-names> Compute median values.\n"
260 " -decilerange=<column-names> Compute decile range.\n"
261 " -percentile=<percentage>,<column-names> Compute specified percentile.\n"
262 " -standarddeviation=<column-names> Compute standard deviations.\n"
263 " -rms=<column-names> Compute RMS values.\n"
264 " -sigma=<column-names> Compute sigma values.\n"
265 " -slope=<indep-parameter>,<column-names> Compute slope for linear fit.\n"
266 " -intercept=<indep-parameter>,<column-names> Compute intercept for linear fit.\n"
267 " -wmean=<weightColumn>,<columnNames> Compute weighted mean.\n"
268 " -wstandarddeviation=<weightColumn>,<columnNames> Compute weighted standard deviation.\n"
269 " -wrms=<weightColumn>,<columnNames> Compute weighted RMS.\n"
270 " -wsigma=<weightColumn>,<columnNames> Compute weighted sigma.\n"
271 " -majorOrder=row|column Set major order.\n\n"
272 "Processes pages from <input> to produce <output> with\n"
273 "one page containing the specified quantities across pages\n"
274 "for each row of the specified columns.\n"
275 "Program by Michael Borland. (" __DATE__ " " __TIME__ ", SVN revision: " SVN_VERSION ")";
276
277int main(int argc, char **argv) {
278 STAT_DEFINITION *stat;
279 long stats;
280 STAT_REQUEST *request;
281 long requests;
282 SCANNED_ARG *scanned; /* structure for scanned arguments */
283 SDDS_DATASET inTable, outTable;
284 long i_arg, code, power, iStat, pages, nowarnings = 0;
285 int64_t i, rows, firstRows;
286 char *input, *output;
287 double *inputData, *otherData, indepData, *weight;
288 unsigned long pipeFlags, majorOrderFlag;
289 double decilePoint[2] = {10.0, 90.0}, decileResult[2];
290 double percentilePoint, percentileResult;
291 double percentile;
292 short columnMajorOrder = -1;
293
295 argc = scanargs(&scanned, argc, argv);
296 if (argc < 2) {
297 bomb("too few arguments", USAGE);
298 exit(EXIT_FAILURE);
299 }
300 weight = NULL;
301 input = output = NULL;
302 stat = NULL;
303 request = NULL;
304 stats = requests = pipeFlags = 0;
305 rows = firstRows = i = 0;
306
307 for (i_arg = 1; i_arg < argc; i_arg++) {
308 if (scanned[i_arg].arg_type == OPTION) {
309 /* process options here */
310 switch (code = match_string(scanned[i_arg].list[0], option, N_OPTIONS, 0)) {
311 case SET_COPY:
312 case SET_MINIMA:
313 case SET_MAXIMA:
314 case SET_LARGEST:
315 case SET_SIGNEDLARGEST:
316 case SET_MEANS:
317 case SET_SDS:
318 case SET_SIGMAS:
319 case SET_RMSS:
320 case SET_MEDIAN:
321 case SET_DRANGE:
322 case SET_EXMM_MEAN:
323 if (scanned[i_arg].n_items < 2) {
324 fprintf(stderr, "error: invalid -%s syntax\n", option[code]);
325 exit(EXIT_FAILURE);
326 }
327 requests = addStatRequests(&request, requests, scanned[i_arg].list + 1, scanned[i_arg].n_items - 1, code, 0, 0, NULL, 0, NULL);
328 break;
329 case SET_WMEANS:
330 case SET_WSDS:
331 case SET_WRMSS:
332 case SET_WSIGMAS:
333 if (scanned[i_arg].n_items < 3) {
334 fprintf(stderr, "error: invalid -%s syntax\n", option[code]);
335 exit(EXIT_FAILURE);
336 }
337 /*note here, items=scanned[i_arg].n_items-2, because the weightedColumn should be excluded */
338 requests = addStatRequests(&request, requests, scanned[i_arg].list + 1, scanned[i_arg].n_items - 1, code, 0, 0, NULL, 1, NULL);
339 break;
340 case SET_SUMS:
341 if (scanned[i_arg].n_items < 3) {
342 fprintf(stderr, "error: invalid -%s syntax\n", option[code]);
343 exit(EXIT_FAILURE);
344 }
345 if (sscanf(scanned[i_arg].list[1], "%ld", &power) != 1 || power < 1) {
346 fprintf(stderr, "error: invalid -%s syntax--bad power in field %s\n", option[code], scanned[i_arg].list[2]);
347 exit(EXIT_FAILURE);
348 }
349 requests = addStatRequests(&request, requests, scanned[i_arg].list + 2, scanned[i_arg].n_items - 2, code, 0, power, NULL, 0, NULL);
350 break;
351 case SET_PERCENTILE:
352 if (scanned[i_arg].n_items < 3) {
353 fprintf(stderr, "error: invalid -%s syntax\n", option[code]);
354 exit(EXIT_FAILURE);
355 }
356 if (sscanf(scanned[i_arg].list[1], "%lf", &percentile) != 1 || percentile < 0 || percentile > 100) {
357 fprintf(stderr, "error: invalid -%s syntax--bad percentage in field %s\n", option[code], scanned[i_arg].list[1]);
358 exit(EXIT_FAILURE);
359 }
360 requests = addStatRequests(&request, requests, scanned[i_arg].list + 2, scanned[i_arg].n_items - 2, code, percentile, 0, NULL, 0, scanned[i_arg].list[1]);
361 break;
362 case SET_SLOPE:
363 case SET_INTERCEPT:
364 case SET_PMINIMA:
365 case SET_PMAXIMA:
366 case SET_CMINIMA:
367 case SET_CMAXIMA:
368 if (scanned[i_arg].n_items < 3) {
369 fprintf(stderr, "error: invalid -%s syntax\n", option[code]);
370 exit(EXIT_FAILURE);
371 }
372 requests = addStatRequests(&request, requests, scanned[i_arg].list + 2, scanned[i_arg].n_items - 2, code, 0, 0, scanned[i_arg].list[1], 0, NULL);
373 break;
374 case SET_PIPE:
375 if (!processPipeOption(scanned[i_arg].list + 1, scanned[i_arg].n_items - 1, &pipeFlags))
376 SDDS_Bomb("invalid -pipe syntax");
377 break;
378 case SET_NOWARNINGS:
379 nowarnings = 1;
380 break;
381 case SET_MAJOR_ORDER:
382 majorOrderFlag = 0;
383 scanned[i_arg].n_items--;
384 if (scanned[i_arg].n_items > 0 && (!scanItemList(&majorOrderFlag, scanned[i_arg].list + 1, &scanned[i_arg].n_items, 0, "row", -1, NULL, 0, SDDS_ROW_MAJOR_ORDER, "column", -1, NULL, 0, SDDS_COLUMN_MAJOR_ORDER, NULL)))
385 SDDS_Bomb("invalid -majorOrder syntax/values");
386 if (majorOrderFlag & SDDS_COLUMN_MAJOR_ORDER)
387 columnMajorOrder = 1;
388 else if (majorOrderFlag & SDDS_ROW_MAJOR_ORDER)
389 columnMajorOrder = 0;
390 break;
391 default:
392 fprintf(stderr, "error: unknown option '%s' given\n", scanned[i_arg].list[0]);
393 exit(EXIT_FAILURE);
394 break;
395 }
396 } else {
397 /* argument is filename */
398 if (!input)
399 input = scanned[i_arg].list[0];
400 else if (!output)
401 output = scanned[i_arg].list[0];
402 else
403 SDDS_Bomb("too many filenames seen");
404 }
405 }
406
407 processFilenames("sddsenvelope", &input, &output, pipeFlags, 0, NULL);
408
409 if (!requests)
410 SDDS_Bomb("no statistics requested");
411
412 if (!SDDS_InitializeInput(&inTable, input))
413 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
414
415 pages = 0;
416 while ((code = SDDS_ReadPage(&inTable)) > 0) {
417 pages++;
418 if (!(rows = SDDS_CountRowsOfInterest(&inTable)))
419 SDDS_Bomb("empty data page in input file");
420 if (code == 1) {
421 firstRows = rows;
422 if (!(stat = compileStatDefinitions(&inTable, &stats, request, requests))) {
423 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors);
424 exit(EXIT_FAILURE);
425 }
426 if (!setupOutputFile(&outTable, output, &inTable, stat, stats, rows, columnMajorOrder)) {
428 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors);
429 else
430 fprintf(stderr, "Error setting up output file.\n");
431 exit(EXIT_FAILURE);
432 }
433 } else if (firstRows != rows)
434 SDDS_Bomb("inconsistent number of rows in input file");
435 for (iStat = 0; iStat < stats; iStat++) {
436 if (stat[iStat].optionCode == SET_COPY) {
437 if (code == 1 && !(stat[iStat].copy = SDDS_GetColumn(&inTable, stat[iStat].sourceColumn)))
438 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
439 continue;
440 }
441 stat[iStat].copy = NULL;
442 if (!(inputData = SDDS_GetColumnInDoubles(&inTable, stat[iStat].sourceColumn)))
443 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
444 switch (stat[iStat].optionCode) {
445 case SET_MINIMA:
446 if (code == 1)
447 for (i = 0; i < rows; i++)
448 stat[iStat].value1[i] = inputData[i];
449 else
450 for (i = 0; i < rows; i++)
451 if (stat[iStat].value1[i] > inputData[i])
452 stat[iStat].value1[i] = inputData[i];
453 break;
454 case SET_MAXIMA:
455 if (code == 1)
456 for (i = 0; i < rows; i++)
457 stat[iStat].value1[i] = inputData[i];
458 else
459 for (i = 0; i < rows; i++)
460 if (stat[iStat].value1[i] < inputData[i])
461 stat[iStat].value1[i] = inputData[i];
462 break;
463 case SET_CMINIMA:
464 /* Value from another column corresponding to minimum value in main column */
465 if (!(otherData = SDDS_GetColumnInDoubles(&inTable, stat[iStat].functionOf)))
466 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
467 if (code == 1)
468 for (i = 0; i < rows; i++) {
469 stat[iStat].value2[i] = inputData[i];
470 stat[iStat].value1[i] = otherData[i];
471 }
472 else
473 for (i = 0; i < rows; i++)
474 if (stat[iStat].value2[i] > inputData[i]) {
475 stat[iStat].value2[i] = inputData[i];
476 stat[iStat].value1[i] = otherData[i];
477 }
478 free(otherData);
479 break;
480 case SET_CMAXIMA:
481 /* Value from another column corresponding to maximum value in main column */
482 if (!(otherData = SDDS_GetColumnInDoubles(&inTable, stat[iStat].functionOf)))
483 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
484 if (code == 1)
485 for (i = 0; i < rows; i++) {
486 stat[iStat].value2[i] = inputData[i];
487 stat[iStat].value1[i] = otherData[i];
488 }
489 else
490 for (i = 0; i < rows; i++)
491 if (stat[iStat].value2[i] < inputData[i]) {
492 stat[iStat].value2[i] = inputData[i];
493 stat[iStat].value1[i] = otherData[i];
494 }
495 free(otherData);
496 break;
497 case SET_LARGEST:
498 if (code == 1)
499 for (i = 0; i < rows; i++)
500 stat[iStat].value1[i] = fabs(inputData[i]);
501 else
502 for (i = 0; i < rows; i++)
503 if (stat[iStat].value1[i] < fabs(inputData[i]))
504 stat[iStat].value1[i] = fabs(inputData[i]);
505 break;
506 case SET_SIGNEDLARGEST:
507 if (code == 1)
508 for (i = 0; i < rows; i++)
509 stat[iStat].value1[i] = inputData[i];
510 else
511 for (i = 0; i < rows; i++)
512 if (fabs(stat[iStat].value1[i]) < fabs(inputData[i]))
513 stat[iStat].value1[i] = inputData[i];
514 break;
515 case SET_MEANS:
516 if (code == 1)
517 for (i = 0; i < rows; i++)
518 stat[iStat].value1[i] = inputData[i];
519 else
520 for (i = 0; i < rows; i++)
521 stat[iStat].value1[i] += inputData[i];
522 break;
523 case SET_WMEANS:
524 if (!(weight = SDDS_GetColumnInDoubles(&inTable, stat[iStat].weightColumn)))
525 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
526 for (i = 0; i < rows; i++) {
527 stat[iStat].sumWeight[i] += weight[i];
528 stat[iStat].value1[i] += inputData[i] * weight[i];
529 }
530 free(weight);
531 break;
532 case SET_SDS:
533 case SET_SIGMAS:
534 if (code == 1)
535 for (i = 0; i < rows; i++) {
536 stat[iStat].value1[i] = inputData[i];
537 stat[iStat].value2[i] = inputData[i] * inputData[i];
538 }
539 else
540 for (i = 0; i < rows; i++) {
541 stat[iStat].value1[i] += inputData[i];
542 stat[iStat].value2[i] += inputData[i] * inputData[i];
543 }
544 break;
545 case SET_WSDS:
546 case SET_WSIGMAS:
547 if (!(weight = SDDS_GetColumnInDoubles(&inTable, stat[iStat].weightColumn)))
548 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
549 for (i = 0; i < rows; i++) {
550 stat[iStat].sumWeight[i] += weight[i];
551 stat[iStat].value1[i] += inputData[i] * weight[i];
552 stat[iStat].value2[i] += inputData[i] * inputData[i] * weight[i];
553 }
554 free(weight);
555 break;
556 case SET_RMSS:
557 if (code == 1)
558 for (i = 0; i < rows; i++)
559 stat[iStat].value1[i] = inputData[i] * inputData[i];
560 else
561 for (i = 0; i < rows; i++)
562 stat[iStat].value1[i] += inputData[i] * inputData[i];
563 break;
564 case SET_WRMSS:
565 if (!(weight = SDDS_GetColumnInDoubles(&inTable, stat[iStat].weightColumn)))
566 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
567 for (i = 0; i < rows; i++) {
568 stat[iStat].sumWeight[i] += weight[i];
569 stat[iStat].value1[i] += inputData[i] * inputData[i] * weight[i];
570 }
571 free(weight);
572 break;
573 case SET_SUMS:
574 if (code == 1)
575 for (i = 0; i < rows; i++)
576 stat[iStat].value1[i] = ipow(inputData[i], stat[iStat].sumPower);
577 else
578 for (i = 0; i < rows; i++)
579 stat[iStat].value1[i] += ipow(inputData[i], stat[iStat].sumPower);
580 break;
581 case SET_PMINIMA:
582 if (!SDDS_GetParameterAsDouble(&inTable, stat[iStat].functionOf, &indepData)) {
583 fprintf(stderr, "error: unable to get value of parameter %s\n", stat[iStat].functionOf);
584 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
585 }
586 if (code == 1)
587 for (i = 0; i < rows; i++) {
588 stat[iStat].value2[i] = inputData[i];
589 stat[iStat].value1[i] = indepData;
590 }
591 else
592 for (i = 0; i < rows; i++) {
593 if (stat[iStat].value2[i] > inputData[i]) {
594 stat[iStat].value2[i] = inputData[i];
595 stat[iStat].value1[i] = indepData;
596 }
597 }
598 break;
599 case SET_PMAXIMA:
600 if (!SDDS_GetParameterAsDouble(&inTable, stat[iStat].functionOf, &indepData)) {
601 fprintf(stderr, "error: unable to get value of parameter %s\n", stat[iStat].functionOf);
602 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
603 }
604 if (code == 1)
605 for (i = 0; i < rows; i++) {
606 stat[iStat].value2[i] = inputData[i];
607 stat[iStat].value1[i] = indepData;
608 }
609 else
610 for (i = 0; i < rows; i++) {
611 if (stat[iStat].value2[i] < inputData[i]) {
612 stat[iStat].value2[i] = inputData[i];
613 stat[iStat].value1[i] = indepData;
614 }
615 }
616 break;
617 case SET_SLOPE:
618 case SET_INTERCEPT:
619 if (!SDDS_GetParameterAsDouble(&inTable, stat[iStat].functionOf, &indepData)) {
620 fprintf(stderr, "error: unable to get value of parameter %s\n", stat[iStat].functionOf);
621 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
622 }
623 /* linear fit:
624 y = a + bx
625 a = (S x^2 Sy - S x S xy)/D
626 b = (N S xy - Sx Sy)/D
627 D = N S x^2 - (S x)^2
628 */
629 if (code == 1)
630 for (i = 0; i < rows; i++) {
631 stat[iStat].value1[i] = indepData; /* Sum x */
632 stat[iStat].value2[i] = indepData * indepData; /* Sum x^2 */
633 stat[iStat].value3[i] = inputData[i]; /* Sum y */
634 stat[iStat].value4[i] = indepData * inputData[i]; /* Sum xy */
635 }
636 else
637 for (i = 0; i < rows; i++) {
638 stat[iStat].value1[i] += indepData;
639 stat[iStat].value2[i] += indepData * indepData;
640 stat[iStat].value3[i] += inputData[i];
641 stat[iStat].value4[i] += indepData * inputData[i];
642 }
643 break;
644 case SET_MEDIAN:
645 case SET_DRANGE:
646 case SET_PERCENTILE:
647 case SET_EXMM_MEAN:
648 if (code == 1)
649 for (i = 0; i < rows; i++) {
650 stat[iStat].array[i] = tmalloc(sizeof(*stat[iStat].array[i]));
651 stat[iStat].array[i][pages - 1] = inputData[i];
652 }
653 else {
654 for (i = 0; i < rows; i++) {
655 stat[iStat].array[i] = SDDS_Realloc(stat[iStat].array[i], sizeof(*stat[iStat].array[i]) * pages);
656 stat[iStat].array[i][pages - 1] = inputData[i];
657 }
658 }
659 break;
660 default:
661 SDDS_Bomb("invalid statistic code (accumulation loop)");
662 break;
663 }
664 free(inputData);
665 }
666 }
667 if (pages == 0)
668 SDDS_Bomb("no pages in input");
669 for (iStat = 0; iStat < stats; iStat++) {
670 switch (stat[iStat].optionCode) {
671 case SET_COPY:
672 case SET_MINIMA:
673 case SET_MAXIMA:
674 case SET_PMINIMA:
675 case SET_PMAXIMA:
676 case SET_CMINIMA:
677 case SET_CMAXIMA:
678 case SET_LARGEST:
679 case SET_SIGNEDLARGEST:
680 case SET_SUMS:
681 break;
682 case SET_MEANS:
683 for (i = 0; i < rows; i++)
684 stat[iStat].value1[i] /= pages;
685 break;
686 case SET_WMEANS:
687 for (i = 0; i < rows; i++)
688 if (stat[iStat].sumWeight[i])
689 stat[iStat].value1[i] /= stat[iStat].sumWeight[i];
690 else {
691 if (!nowarnings)
692 fprintf(stderr, "warning: the total weight for the %" PRId64 "th row of %s is zero.\n", i + 1, stat[iStat].sourceColumn);
693 stat[iStat].value1[i] = DBL_MAX;
694 }
695 break;
696 case SET_SDS:
697 if (pages < 2)
698 stat[iStat].value1[i] = DBL_MAX;
699 else
700 for (i = 0; i < rows; i++) {
701 double tmp1;
702 if ((tmp1 = stat[iStat].value2[i] / pages - sqr(stat[iStat].value1[i] / pages)) <= 0)
703 stat[iStat].value1[i] = 0;
704 else
705 stat[iStat].value1[i] = sqrt(tmp1 * pages / (pages - 1.0));
706 }
707 break;
708 case SET_WSDS:
709 if (pages < 2)
710 stat[iStat].value1[i] = DBL_MAX;
711 else
712 for (i = 0; i < rows; i++) {
713 double tmp1;
714 if (stat[iStat].sumWeight[i]) {
715 if ((tmp1 = stat[iStat].value2[i] / stat[iStat].sumWeight[i] - sqr(stat[iStat].value1[i] / stat[iStat].sumWeight[i])) <= 0)
716 stat[iStat].value1[i] = 0;
717 else
718 stat[iStat].value1[i] = sqrt(tmp1 * pages / (pages - 1.0));
719 } else {
720 if (!nowarnings)
721 fprintf(stderr, "Warning, the total weight for the %" PRId64 "th row of %s is zero.\n", i + 1, stat[iStat].sourceColumn);
722 stat[iStat].value1[i] = DBL_MAX;
723 }
724 }
725 break;
726 case SET_SIGMAS:
727 if (pages < 2)
728 stat[iStat].value1[i] = DBL_MAX;
729 else
730 for (i = 0; i < rows; i++) {
731 double tmp1;
732 if ((tmp1 = stat[iStat].value2[i] / pages - sqr(stat[iStat].value1[i] / pages)) <= 0)
733 stat[iStat].value1[i] = 0;
734 else
735 stat[iStat].value1[i] = sqrt(tmp1 / (pages - 1.0));
736 }
737 break;
738 case SET_WSIGMAS:
739 if (pages < 2)
740 stat[iStat].value1[i] = DBL_MAX;
741 else
742 for (i = 0; i < rows; i++) {
743 double tmp1;
744 if (stat[iStat].sumWeight[i]) {
745 if ((tmp1 = stat[iStat].value2[i] / stat[iStat].sumWeight[i] - sqr(stat[iStat].value1[i] / stat[iStat].sumWeight[i])) <= 0)
746 stat[iStat].value1[i] = 0;
747 else
748 stat[iStat].value1[i] = sqrt(tmp1 / (pages - 1.0));
749 } else {
750 if (!nowarnings)
751 fprintf(stderr, "Warning, the total weight for the %" PRId64 "th row of %s is zero.\n", i + 1, stat[iStat].sourceColumn);
752 stat[iStat].value1[i] = DBL_MAX;
753 }
754 }
755 break;
756 case SET_RMSS:
757 for (i = 0; i < rows; i++)
758 stat[iStat].value1[i] = sqrt(stat[iStat].value1[i] / pages);
759 break;
760 case SET_WRMSS:
761 for (i = 0; i < rows; i++) {
762 if (stat[iStat].sumWeight[i])
763 stat[iStat].value1[i] = sqrt(stat[iStat].value1[i] / stat[iStat].sumWeight[i]);
764 else {
765 if (!nowarnings)
766 fprintf(stderr, "Warning, the total weight for the %" PRId64 "th row of %s is zero.\n", i + 1, stat[iStat].sourceColumn);
767 stat[iStat].value1[i] = DBL_MAX;
768 }
769 }
770 break;
771 case SET_SLOPE:
772 for (i = 0; i < rows; i++) {
773 double D;
774 D = pages * stat[iStat].value2[i] - stat[iStat].value1[i] * stat[iStat].value1[i];
775 stat[iStat].value1[i] = (pages * stat[iStat].value4[i] - stat[iStat].value1[i] * stat[iStat].value3[i]) / D;
776 }
777 break;
778 case SET_INTERCEPT:
779 for (i = 0; i < rows; i++) {
780 double D;
781 D = pages * stat[iStat].value2[i] - stat[iStat].value1[i] * stat[iStat].value1[i];
782 stat[iStat].value1[i] = (stat[iStat].value2[i] * stat[iStat].value3[i] - stat[iStat].value1[i] * stat[iStat].value4[i]) / D;
783 }
784 break;
785 case SET_MEDIAN:
786 for (i = 0; i < rows; i++) {
787 compute_median(&stat[iStat].value1[i], stat[iStat].array[i], pages);
788 }
789 break;
790 case SET_DRANGE:
791 for (i = 0; i < rows; i++) {
792 if (!compute_percentiles(decileResult, decilePoint, 2, stat[iStat].array[i], pages))
793 stat[iStat].value1[i] = 0;
794 else
795 stat[iStat].value1[i] = decileResult[1] - decileResult[0];
796 }
797 break;
798 case SET_PERCENTILE:
799 percentilePoint = stat[iStat].percentile;
800 for (i = 0; i < rows; i++) {
801 if (!compute_percentiles(&percentileResult, &percentilePoint, 1, stat[iStat].array[i], pages))
802 stat[iStat].value1[i] = 0;
803 else
804 stat[iStat].value1[i] = percentileResult;
805 }
806 break;
807 case SET_EXMM_MEAN:
808 for (i = 0; i < rows; i++)
809 if (!compute_mean_exclude_min_max(&(stat[iStat].value1[i]), stat[iStat].array[i], pages))
810 stat[iStat].value1[i] = 0;
811 break;
812 default:
813 SDDS_Bomb("invalid statistic code (final loop)");
814 break;
815 }
816 if (stat[iStat].optionCode == SET_COPY) {
817 if (!SDDS_SetColumn(&outTable, SDDS_SET_BY_NAME, stat[iStat].copy, rows, stat[iStat].resultColumn)) {
818 fprintf(stderr, "error setting column values for column %s\n", stat[iStat].resultColumn);
819 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
820 }
821 } else if (!SDDS_SetColumnFromDoubles(&outTable, SDDS_SET_BY_NAME, stat[iStat].value1, rows, stat[iStat].resultColumn)) {
822 fprintf(stderr, "error setting column values for column %s\n", stat[iStat].resultColumn);
823 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
824 }
825 if (stat[iStat].value1)
826 free(stat[iStat].value1);
827 if (stat[iStat].value2)
828 free(stat[iStat].value2);
829 if (stat[iStat].value3)
830 free(stat[iStat].value3);
831 if (stat[iStat].value4)
832 free(stat[iStat].value4);
833 if (stat[iStat].copy)
834 free(stat[iStat].copy);
835 if (stat[iStat].array) {
836 for (i = 0; i < rows; i++) {
837 free(stat[iStat].array[i]);
838 }
839 free(stat[iStat].array);
840 }
841 if (stat[iStat].sumWeight)
842 free(stat[iStat].sumWeight);
843 free(stat[iStat].sourceColumn);
844 free(stat[iStat].resultColumn);
845 stat[iStat].value1 = stat[iStat].value2 = stat[iStat].value3 = stat[iStat].value4 = NULL;
846 stat[iStat].copy = NULL;
847 stat[iStat].array = NULL;
848 }
849 free(stat);
850 if (!SDDS_WritePage(&outTable) || !SDDS_Terminate(&inTable) || !SDDS_Terminate(&outTable))
851 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
852 return EXIT_SUCCESS;
853}
854
855long addStatRequests(STAT_REQUEST **statRequest, long requests, char **item, long items, long code, double percentile, long power, char *functionOf, long weighted, char *percentileString) {
856 long i;
857 /*weighted factor should be either 0 or 1 */
858 if (weighted != 0 && weighted != 1)
859 SDDS_Bomb("addStatRequests: weighted parameter should be either 0 or 1");
860 if (code == SET_PERCENTILE && (!percentileString || !strlen(percentileString))) {
861 fprintf(stderr, "Percentile specification is incorrect: percentile=%e, percentileString=%s\n", percentile, percentileString ? percentileString : "NULL");
862 exit(EXIT_FAILURE);
863 }
864 *statRequest = SDDS_Realloc(*statRequest, sizeof(**statRequest) * (requests + items - weighted));
865 for (i = 0; i < items - weighted; i++) {
866 if (weighted)
867 (*statRequest)[requests + i].weightColumnName = item[0];
868 else
869 (*statRequest)[requests + i].weightColumnName = NULL;
870 (*statRequest)[i + requests].columnName = item[i + weighted];
871 (*statRequest)[i + requests].optionCode = code;
872 (*statRequest)[i + requests].sumPower = power;
873 (*statRequest)[i + requests].percentile = percentile;
874 (*statRequest)[i + requests].percentileString = percentileString;
875 (*statRequest)[i + requests].functionOf = functionOf;
876 }
877
878 return items + requests - weighted;
879}
880
881STAT_DEFINITION *compileStatDefinitions(SDDS_DATASET *inTable, long *stats, STAT_REQUEST *request, long requests) {
882 STAT_DEFINITION *stat;
883 long iReq, iStat, iName;
884 int32_t columnNames;
885 char s[SDDS_MAXLINE];
886 char **columnName;
887
888 *stats = iStat = 0;
889 stat = tmalloc(sizeof(*stat) * requests);
890 for (iReq = 0; iReq < requests; iReq++) {
891 if (iStat >= *stats)
892 stat = SDDS_Realloc(stat, sizeof(*stat) * (*stats += 10));
893 if (!has_wildcards(request[iReq].columnName)) {
894 if (SDDS_GetColumnIndex(inTable, request[iReq].columnName) < 0) {
895 sprintf(s, "error: column %s not found input file", request[iReq].columnName);
896 SDDS_SetError(s);
897 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
898 }
899 stat[iStat].weightColumn = request[iReq].weightColumnName;
900 stat[iStat].sourceColumn = request[iReq].columnName;
901 stat[iStat].optionCode = request[iReq].optionCode;
902 stat[iStat].percentile = request[iReq].percentile;
903 stat[iStat].percentileString = request[iReq].percentileString;
904 stat[iStat].sumPower = request[iReq].sumPower;
905 stat[iStat].value1 = stat[iStat].value2 = stat[iStat].value3 = stat[iStat].value4 = NULL;
906 stat[iStat].array = NULL;
907 stat[iStat].copy = NULL;
908 stat[iStat].sumWeight = NULL;
909 if ((stat[iStat].functionOf = request[iReq].functionOf)) {
910 if (stat[iStat].optionCode != SET_CMAXIMA && stat[iStat].optionCode != SET_CMINIMA) {
911 if (SDDS_GetParameterIndex(inTable, request[iReq].functionOf) < 0) {
912 sprintf(s, "error: parameter %s not found input file (1)", request[iReq].functionOf);
913 SDDS_SetError(s);
914 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
915 }
916 } else {
917 if (SDDS_GetColumnIndex(inTable, request[iReq].functionOf) < 0) {
918 sprintf(s, "error: column %s not found input file (1)", request[iReq].functionOf);
919 SDDS_SetError(s);
920 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
921 }
922 }
923 }
924 iStat++;
925 } else {
926 SDDS_SetColumnFlags(inTable, 0);
927 if (!SDDS_SetColumnsOfInterest(inTable, SDDS_MATCH_STRING, request[iReq].columnName, SDDS_OR))
928 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
929 if (!(columnName = SDDS_GetColumnNames(inTable, &columnNames))) {
930 sprintf(s, "no columns selected for wildcard sequence %s", request[iReq].columnName);
931 SDDS_SetError(s);
932 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
933 }
934 if (iStat + columnNames > *stats)
935 stat = SDDS_Realloc(stat, sizeof(*stat) * (*stats = iStat + columnNames + 10));
936 for (iName = 0; iName < columnNames; iName++) {
937 stat[iStat + iName].weightColumn = request[iReq].weightColumnName;
938 stat[iStat + iName].sourceColumn = columnName[iName];
939 stat[iStat + iName].optionCode = request[iReq].optionCode;
940 stat[iStat + iName].sumPower = request[iReq].sumPower;
941 stat[iStat + iName].percentile = request[iReq].percentile;
942 stat[iStat + iName].percentileString = request[iReq].percentileString;
943 stat[iStat + iName].value1 = stat[iStat + iName].value2 = stat[iStat + iName].value3 = stat[iStat + iName].value4 = NULL;
944 stat[iStat + iName].array = NULL;
945 stat[iStat + iName].copy = NULL;
946 stat[iStat + iName].sumWeight = NULL;
947 if ((stat[iStat + iName].functionOf = request[iReq].functionOf) && iName == 0) {
948 if (stat[iStat + iName].optionCode != SET_CMAXIMA && stat[iStat + iName].optionCode != SET_CMINIMA) {
949 if (SDDS_GetParameterIndex(inTable, request[iReq].functionOf) < 0) {
950 sprintf(s, "error: parameter %s not found input file (2)", request[iReq].functionOf);
951 SDDS_SetError(s);
952 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
953 }
954 } else {
955 if (SDDS_GetColumnIndex(inTable, request[iReq].functionOf) < 0) {
956 sprintf(s, "error: column %s not found input file (2)", request[iReq].functionOf);
957 SDDS_SetError(s);
958 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
959 }
960 }
961 }
962 }
963 iStat += columnNames;
964 free(columnName);
965 }
966 }
967
968 *stats = iStat;
969 for (iStat = 0; iStat < *stats; iStat++) {
970 switch (stat[iStat].optionCode) {
971 case SET_COPY:
972 strcpy(s, stat[iStat].sourceColumn);
973 break;
974 case SET_SUMS:
975 if (stat[iStat].sumPower == 1)
976 sprintf(s, "%s%s", stat[iStat].sourceColumn, optionSuffix[stat[iStat].optionCode]);
977 else
978 sprintf(s, "%s%ld%s", stat[iStat].sourceColumn, stat[iStat].sumPower, optionSuffix[stat[iStat].optionCode]);
979 break;
980 case SET_PERCENTILE:
981 sprintf(s, "%s%s%s", stat[iStat].sourceColumn, stat[iStat].percentileString, optionSuffix[stat[iStat].optionCode]);
982 break;
983 case SET_PMAXIMA:
984 case SET_PMINIMA:
985 sprintf(s, "%s%s%s", stat[iStat].functionOf, optionSuffix[stat[iStat].optionCode], stat[iStat].sourceColumn);
986 break;
987 case SET_CMAXIMA:
988 case SET_CMINIMA:
989 sprintf(s, "%s%s%s", stat[iStat].functionOf, optionSuffix[stat[iStat].optionCode], stat[iStat].sourceColumn);
990 break;
991 default:
992 sprintf(s, "%s%s", stat[iStat].sourceColumn, optionSuffix[stat[iStat].optionCode]);
993 break;
994 }
995 if (!SDDS_CopyString(&stat[iStat].resultColumn, s))
996 SDDS_PrintErrors(stderr, SDDS_VERBOSE_PrintErrors | SDDS_EXIT_PrintErrors);
997 }
998 return stat;
999}
1000
1001long setupOutputFile(SDDS_DATASET *outTable, char *output, SDDS_DATASET *inTable, STAT_DEFINITION *stat, long stats, int64_t rows, short columnMajorOrder) {
1002 long column;
1003 char s[SDDS_MAXLINE], *symbol, *symbol1, *units1;
1004
1005 if (!SDDS_InitializeOutput(outTable, SDDS_BINARY, 0, NULL, "sddsenvelope output", output))
1006 return 0;
1007 if (columnMajorOrder != -1)
1008 outTable->layout.data_mode.column_major = columnMajorOrder;
1009 else
1010 outTable->layout.data_mode.column_major = inTable->layout.data_mode.column_major;
1011 for (column = 0; column < stats; column++) {
1012 stat[column].value1 = calloc(sizeof(*stat[column].value1), rows);
1013 stat[column].value2 = stat[column].value3 = stat[column].value4 = NULL;
1014 if (stat[column].optionCode == SET_SDS || stat[column].optionCode == SET_SIGMAS ||
1015 stat[column].optionCode == SET_WSDS || stat[column].optionCode == SET_WSIGMAS ||
1016 stat[column].optionCode == SET_PMINIMA || stat[column].optionCode == SET_PMAXIMA ||
1017 stat[column].optionCode == SET_CMINIMA || stat[column].optionCode == SET_CMAXIMA)
1018 stat[column].value2 = calloc(sizeof(*stat[column].value2), rows);
1019 if (stat[column].optionCode == SET_INTERCEPT || stat[column].optionCode == SET_SLOPE) {
1020 stat[column].value2 = malloc(sizeof(*stat[column].value2) * rows);
1021 stat[column].value3 = malloc(sizeof(*stat[column].value3) * rows);
1022 stat[column].value4 = malloc(sizeof(*stat[column].value4) * rows);
1023 }
1024 if (stat[column].optionCode == SET_WSDS || stat[column].optionCode == SET_WSIGMAS ||
1025 stat[column].optionCode == SET_WRMSS || stat[column].optionCode == SET_WMEANS)
1026 stat[column].sumWeight = calloc(sizeof(*stat[column].sumWeight), rows);
1027 if (stat[column].optionCode == SET_MEDIAN || stat[column].optionCode == SET_DRANGE ||
1028 stat[column].optionCode == SET_PERCENTILE || stat[column].optionCode == SET_EXMM_MEAN) {
1029 stat[column].array = tmalloc(sizeof(*stat[column].array) * rows);
1030 }
1031 if (!SDDS_TransferColumnDefinition(outTable, inTable, stat[column].sourceColumn, stat[column].resultColumn)) {
1032 sprintf(s, "Problem transferring definition of column %s to %s\n", stat[column].sourceColumn, stat[column].resultColumn);
1033 SDDS_SetError(s);
1034 return 0;
1035 }
1036 if (SDDS_ChangeColumnInformation(outTable, "description", NULL, SDDS_SET_BY_NAME, stat[column].resultColumn) != SDDS_STRING ||
1037 SDDS_GetColumnInformation(outTable, "symbol", &symbol, SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1038 fprintf(stderr, "Error: problem setting description for column %s\n", stat[column].resultColumn);
1039 return 0;
1040 }
1041 if (stat[column].optionCode > 0) {
1042 if (SDDS_ChangeColumnInformation(outTable, "type", "double", SDDS_PASS_BY_STRING | SDDS_SET_BY_NAME, stat[column].resultColumn) != SDDS_LONG) {
1043 fprintf(stderr, "Error: problem setting type for column %s\n", stat[column].resultColumn);
1044 return 0;
1045 }
1046 }
1047 if (!symbol)
1048 SDDS_CopyString(&symbol, stat[column].sourceColumn);
1049 switch (stat[column].optionCode) {
1050 case SET_COPY:
1051 strcpy(s, symbol);
1052 break;
1053 case SET_SUMS:
1054 if (stat[column].sumPower == 1)
1055 sprintf(s, "%s[%s]", optionSuffix[stat[column].optionCode], symbol);
1056 else
1057 sprintf(s, "%s[%s$a%ld$n]", optionSuffix[stat[column].optionCode], symbol, stat[column].sumPower);
1058 break;
1059 case SET_PERCENTILE:
1060 sprintf(s, "%s[%s,%g]", optionSuffix[stat[column].optionCode], symbol, stat[column].percentile);
1061 break;
1062 case SET_PMINIMA:
1063 case SET_PMAXIMA:
1064 if (SDDS_GetParameterInformation(inTable, "symbol", &symbol1, SDDS_BY_NAME, stat[column].functionOf) != SDDS_STRING ||
1065 !symbol1 ||
1066 !strlen(symbol1))
1067 symbol1 = stat[column].functionOf;
1068 sprintf(s, "%s[%s:%s]", optionSuffix[stat[column].optionCode], symbol, symbol1);
1069 if (SDDS_GetParameterInformation(inTable, "units", &units1, SDDS_BY_NAME, stat[column].functionOf) != SDDS_STRING)
1070 return 0;
1071 if (units1) {
1072 if (SDDS_ChangeColumnInformation(outTable, "units", units1, SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1073 fprintf(stderr, "Error: problem setting units for column %s (1)\n", stat[column].resultColumn);
1074 return 0;
1075 }
1076 } else if (SDDS_ChangeColumnInformation(outTable, "units", "", SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1077 fprintf(stderr, "Error: problem setting units for column %s (2)\n", stat[column].resultColumn);
1078 return 0;
1079 }
1080 break;
1081 case SET_CMINIMA:
1082 case SET_CMAXIMA:
1083 if (SDDS_GetColumnInformation(inTable, "symbol", &symbol1, SDDS_BY_NAME, stat[column].functionOf) != SDDS_STRING ||
1084 !symbol1 ||
1085 !strlen(symbol1))
1086 symbol1 = stat[column].functionOf;
1087 sprintf(s, "%s[%s:%s]", optionSuffix[stat[column].optionCode], symbol, symbol1);
1088 if (SDDS_GetColumnInformation(inTable, "units", &units1, SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING)
1089 return 0;
1090 if (units1) {
1091 if (SDDS_ChangeColumnInformation(outTable, "units", units1, SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1092 fprintf(stderr, "Error: problem setting units for column %s\n", stat[column].resultColumn);
1093 return 0;
1094 }
1095 } else if (SDDS_ChangeColumnInformation(outTable, "units", "", SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1096 fprintf(stderr, "Error: problem setting units for column %s\n", stat[column].resultColumn);
1097 return 0;
1098 }
1099 break;
1100 case SET_INTERCEPT:
1101 case SET_SLOPE:
1102 if (SDDS_GetParameterInformation(inTable, "symbol", &symbol1, SDDS_BY_NAME, stat[column].functionOf) != SDDS_STRING ||
1103 !symbol1 ||
1104 !strlen(symbol1))
1105 symbol1 = stat[column].functionOf;
1106 sprintf(s, "%s[%s:%s]", optionSuffix[stat[column].optionCode], symbol, symbol1);
1107 break;
1108 default:
1109 sprintf(s, "%s[%s]", optionSuffix[stat[column].optionCode], symbol);
1110 break;
1111 }
1112 free(symbol);
1113 if (SDDS_ChangeColumnInformation(outTable, "symbol", s, SDDS_BY_NAME, stat[column].resultColumn) != SDDS_STRING) {
1114 fprintf(stderr, "Error: problem setting symbol for column %s\n", stat[column].resultColumn);
1115 return 0;
1116 }
1117 }
1118 if (!SDDS_WriteLayout(outTable) || !SDDS_StartPage(outTable, rows))
1119 return 0;
1120 return 1;
1121}
1122
1123int compute_mean_exclude_min_max(double *value, double *data, long n) {
1124 double sum = 0;
1125 long i, count = 0;
1126 double min, max;
1127 max = -(min = DBL_MAX);
1128 if (n <= 0 || !data || !value)
1129 return 0;
1130 for (i = 0; i < n; i++) {
1131 if (min > data[i])
1132 min = data[i];
1133 if (max < data[i])
1134 max = data[i];
1135 }
1136 for (i = 0; i < n; i++) {
1137 if (data[i] == min || data[i] == max)
1138 continue;
1139 count++;
1140 sum += data[i];
1141 }
1142 if (count == 0)
1143 *value = min;
1144 else
1145 *value = sum / count;
1146 return 1;
1147}
SDDS (Self Describing Data Set) Data Types Definitions and Function Prototypes.
int32_t SDDS_StartPage(SDDS_DATASET *SDDS_dataset, int64_t expected_n_rows)
int32_t SDDS_SetColumnFromDoubles(SDDS_DATASET *SDDS_dataset, int32_t mode, double *data, int64_t rows,...)
Sets the values for a single data column using double-precision floating-point numbers.
int32_t SDDS_SetColumn(SDDS_DATASET *SDDS_dataset, int32_t mode, void *data, int64_t rows,...)
Sets the values for one data column in the current data table of an SDDS dataset.
void * SDDS_GetColumn(SDDS_DATASET *SDDS_dataset, char *column_name)
Retrieves a copy of the data for a specified column, including only rows marked as "of interest".
double * SDDS_GetParameterAsDouble(SDDS_DATASET *SDDS_dataset, char *parameter_name, double *memory)
Retrieves the value of a specified parameter as a double from the current data table of an SDDS datas...
int64_t SDDS_CountRowsOfInterest(SDDS_DATASET *SDDS_dataset)
Counts the number of rows marked as "of interest" in the current data table.
int32_t SDDS_SetColumnsOfInterest(SDDS_DATASET *SDDS_dataset, int32_t mode,...)
Sets the acceptance flags for columns based on specified naming criteria.
int32_t SDDS_SetColumnFlags(SDDS_DATASET *SDDS_dataset, int32_t column_flag_value)
Sets the acceptance flags for all columns in the current data table of a data set.
double * SDDS_GetColumnInDoubles(SDDS_DATASET *SDDS_dataset, char *column_name)
Retrieves the data of a specified numerical column as an array of doubles, considering only rows mark...
int32_t SDDS_GetParameterInformation(SDDS_DATASET *SDDS_dataset, char *field_name, void *memory, int32_t mode,...)
Retrieves information about a specified parameter in the SDDS dataset.
Definition SDDS_info.c:117
int32_t SDDS_ChangeColumnInformation(SDDS_DATASET *SDDS_dataset, char *field_name, void *memory, int32_t mode,...)
Modifies a specific field in a column definition within the SDDS dataset.
Definition SDDS_info.c:364
int32_t SDDS_GetColumnInformation(SDDS_DATASET *SDDS_dataset, char *field_name, void *memory, int32_t mode,...)
Retrieves information about a specified column in the SDDS dataset.
Definition SDDS_info.c:41
int32_t SDDS_InitializeInput(SDDS_DATASET *SDDS_dataset, char *filename)
Definition SDDS_input.c:49
int32_t SDDS_Terminate(SDDS_DATASET *SDDS_dataset)
int32_t SDDS_ReadPage(SDDS_DATASET *SDDS_dataset)
int32_t SDDS_InitializeOutput(SDDS_DATASET *SDDS_dataset, int32_t data_mode, int32_t lines_per_row, const char *description, const char *contents, const char *filename)
Initializes the SDDS output dataset.
int32_t SDDS_WritePage(SDDS_DATASET *SDDS_dataset)
Writes the current data table to the output file.
int32_t SDDS_WriteLayout(SDDS_DATASET *SDDS_dataset)
Writes the SDDS layout header to the output file.
int32_t SDDS_TransferColumnDefinition(SDDS_DATASET *target, SDDS_DATASET *source, char *name, char *newName)
Transfers a column definition from a source dataset to a target dataset.
void SDDS_SetError(char *error_text)
Records an error message in the SDDS error stack.
Definition SDDS_utils.c:379
int32_t SDDS_GetParameterIndex(SDDS_DATASET *SDDS_dataset, char *name)
Retrieves the index of a named parameter in the SDDS dataset.
int32_t SDDS_GetColumnIndex(SDDS_DATASET *SDDS_dataset, char *name)
Retrieves the index of a named column in the SDDS dataset.
char ** SDDS_GetColumnNames(SDDS_DATASET *SDDS_dataset, int32_t *number)
Retrieves the names of all columns in the SDDS dataset.
void SDDS_PrintErrors(FILE *fp, int32_t mode)
Prints recorded error messages to a specified file stream.
Definition SDDS_utils.c:432
void SDDS_RegisterProgramName(const char *name)
Registers the executable program name for use in error messages.
Definition SDDS_utils.c:288
int32_t SDDS_NumberOfErrors()
Retrieves the number of errors recorded by SDDS library routines.
Definition SDDS_utils.c:304
void SDDS_Bomb(char *message)
Terminates the program after printing an error message and recorded errors.
Definition SDDS_utils.c:342
int32_t SDDS_CopyString(char **target, const char *source)
Copies a source string to a target string with memory allocation.
Definition SDDS_utils.c:856
void * SDDS_Realloc(void *old_ptr, size_t new_size)
Reallocates memory to a new size.
Definition SDDS_utils.c:677
#define SDDS_STRING
Identifier for the string data type.
Definition SDDStypes.h:85
#define SDDS_LONG
Identifier for the signed 32-bit integer data type.
Definition SDDStypes.h:61
void * tmalloc(uint64_t size_of_block)
Allocates a memory block of the specified size with zero initialization.
Definition array.c:59
void bomb(char *error, char *usage)
Reports error messages to the terminal and aborts the program.
Definition bomb.c:26
double ipow(const double x, const int64_t p)
Compute x raised to the power p (x^p).
Definition ipow.c:33
long match_string(char *string, char **option, long n_options, long mode)
Matches a given string against an array of option strings based on specified modes.
long compute_percentiles(double *position, double *percent, long positions, double *x, long n)
Computes multiple percentiles of an array of doubles.
Definition median.c:84
long compute_median(double *value, double *x, long n)
Computes the median of an array of doubles.
Definition median.c:29
int scanargs(SCANNED_ARG **scanned, int argc, char **argv)
Definition scanargs.c:36
long processPipeOption(char **item, long items, unsigned long *flags)
Definition scanargs.c:356
void processFilenames(char *programName, char **input, char **output, unsigned long pipeFlags, long noWarnings, long *tmpOutputUsed)
Definition scanargs.c:390
long scanItemList(unsigned long *flags, char **item, long *items, unsigned long mode,...)
Scans a list of items and assigns values based on provided keywords and types.
int has_wildcards(char *template)
Check if a template string contains any wildcard characters.
Definition wild_match.c:498