From 09f7f3ec8bbff5102d0447ee664df3f3a404c9bc Mon Sep 17 00:00:00 2001 From: tlatorre Date: Fri, 12 Oct 2018 10:20:28 -0500 Subject: skip PMTs which weren't hit for the fast likelihood calculation --- src/fit.c | 8 ++++++-- src/likelihood.c | 18 +++++++++++++++++- 2 files changed, 23 insertions(+), 3 deletions(-) diff --git a/src/fit.c b/src/fit.c index 4d8d389..3ab0ffe 100644 --- a/src/fit.c +++ b/src/fit.c @@ -5253,8 +5253,6 @@ int fit_event(event *ev, double *xopt, double *fmin) } } - memcpy(x,xopt,sizeof(x)); - /* Reset the lower and upper bounds. */ lb[0] = Tmin; lb[1] = -1000.0; @@ -5286,6 +5284,12 @@ int fit_event(event *ev, double *xopt, double *fmin) nlopt_set_ftol_abs(opt, 1e-5); nlopt_set_maxeval(opt, 1000); + memcpy(x,xopt,sizeof(x)); + + nlopt_optimize(opt,x,&fval); + + if (stop) goto end; + do { *fmin = fval; memcpy(xopt,x,sizeof(x)); diff --git a/src/likelihood.c b/src/likelihood.c index ccdfc84..5f66a72 100644 --- a/src/likelihood.c +++ b/src/likelihood.c @@ -392,6 +392,10 @@ double nll_muon(event *ev, double T0, double *pos, double *dir, double t0, doubl for (i = 0; i < MAX_PMTS; i++) { if (pmts[i].pmt_type != PMT_NORMAL) continue; + /* Skip PMTs which weren't hit when doing the "fast" likelihood + * calculation. */ + if (fast && !ev->pmt_hits[i].hit) continue; + params.i = i; if (fast) { @@ -501,13 +505,25 @@ double nll_muon(event *ev, double T0, double *pos, double *dir, double t0, doubl for (i = 0; i < MAX_PMTS; i++) { if (ev->pmt_hits[i].flags || pmts[i].pmt_type != PMT_NORMAL) continue; - mu[i] = mu_direct[i] + mu_indirect + mu_noise; + + /* Skip PMTs which weren't hit when doing the "fast" likelihood + * calculation. */ + if (fast && !ev->pmt_hits[i].hit) continue; + + if (fast) + mu[i] = mu_direct[i] + mu_noise; + else + mu[i] = mu_direct[i] + mu_indirect + mu_noise; } nhit = 0; for (i = 0; i < MAX_PMTS; i++) { if (ev->pmt_hits[i].flags || pmts[i].pmt_type != PMT_NORMAL) continue; + /* Skip PMTs which weren't hit when doing the "fast" likelihood + * calculation. */ + if (fast && !ev->pmt_hits[i].hit) continue; + log_mu = log(mu[i]); if (ev->pmt_hits[i].hit) { -- cgit