diff options
Diffstat (limited to 'src/fit.c')
-rw-r--r-- | src/fit.c | 73 |
1 files changed, 66 insertions, 7 deletions
@@ -12,6 +12,9 @@ #include <nlopt.h> #include <math.h> /* for sin(), cos(), etc. */ #include <sys/time.h> /* for gettimeofday() */ +#include <inttypes.h> /* for PRIu32 macro */ +#include <string.h> /* for memcpy() */ +#include <errno.h> /* for errno */ #define EV_RECORD 0x45562020 // 'EV ' (as written to ZDAB file) @@ -64,7 +67,7 @@ double nll(unsigned int n, const double *x, double *grad, void *params) return fval; } -int fit_event(event *ev, double *fmin) +int fit_event(event *ev, double *xopt, double *fmin) { double x[9], ss[9], lb[9], ub[9]; int rv; @@ -120,11 +123,21 @@ int fit_event(event *ev, double *fmin) rv = nlopt_optimize(opt,x,fmin); + memcpy(xopt,x,sizeof(x)); + nlopt_destroy(opt); return rv; } +void usage(void) +{ + fprintf(stderr,"Usage: ./fit [options] FILENAME\n"); + fprintf(stderr," -o output file\n"); + fprintf(stderr," -h display this help message\n"); + exit(1); +} + int main(int argc, char **argv) { int i; @@ -137,19 +150,47 @@ int main(int argc, char **argv) int crate, card, channel; int id; double fmin; - - if (argc < 2) { - fprintf(stderr, "usage: fit [filename]\n"); - return 1; + double xopt[9]; + char *filename = NULL; + char *output = NULL; + FILE *fout = NULL; + + for (i = 1; i < argc; i++) { + if (argv[i][0] == '-') { + switch (argv[i][1]) { + case 'o': + output = argv[++i]; + break; + case 'h': + usage(); + default: + fprintf(stderr, "unrecognized option '%s'\n", argv[i]); + exit(1); + } + } else { + filename = argv[i]; + } } - f = zebra_open(argv[1]); + if (!filename) + usage(); + + f = zebra_open(filename); if (!f) { fprintf(stderr, "%s\n", zebra_err); return 1; } + if (output) { + fout = fopen(output, "w"); + + if (!fout) { + fprintf(stderr, "failed to open '%s': %s\n", output, strerror(errno)); + return 1; + } + } + load_pmt_info(); for (i = 0; i < MAX_PMTS; i++) { @@ -201,7 +242,23 @@ int main(int argc, char **argv) if (ev.run != -1) { if (ev.run != bev.run || ev.gtid != bev.gtr_id) { /* New event, so we need to fit the old event. */ - fit_event(&ev,&fmin); + fit_event(&ev,xopt,&fmin); + + if (fout) { + fprintf(fout, "%" PRIu32 " %10.2f %7.2f %7.2f %7.2f %4.2f %4.2f %5.2f %5.2f %5.2f f() = %7.3e\n", + ev.gtid, + xopt[0], + xopt[1], + xopt[2], + xopt[3], + xopt[4], + xopt[5], + xopt[6], + xopt[7], + xopt[8], + fmin); + fflush(fout); + } } } ev.run = bev.run; @@ -217,6 +274,8 @@ int main(int argc, char **argv) db_free(db); + if (fout) fclose(fout); + zebra_close(f); return 0; |