diff options
Diffstat (limited to 'src/fit.c')
| -rw-r--r-- | src/fit.c | 20 |
1 files changed, 13 insertions, 7 deletions
@@ -5298,11 +5298,6 @@ int bisect_energy(double qhs_sum, double distance, int id, double cos_theta, dou pars.qhs_sum = qhs_sum; pars.cos_theta = cos_theta; - s = gsl_root_fsolver_alloc(gsl_root_fsolver_brent); - - F.function = &get_expected_photons_root; - F.params = &pars; - switch (id) { case IDP_E_MINUS: case IDP_E_PLUS: @@ -5326,6 +5321,17 @@ int bisect_energy(double qhs_sum, double distance, int id, double cos_theta, dou /* Calculate the Cerenkov threshold for a muon. */ Tmin = mass/sqrt(1.0-1.0/(n_d2o*n_d2o)) - mass; + if (get_expected_photons(id,Tmin,distance,cos_theta) > qhs_sum) + return Tmin; + + if (get_expected_photons(id,Tmax,distance,cos_theta) < qhs_sum) + return Tmax; + + s = gsl_root_fsolver_alloc(gsl_root_fsolver_brent); + + F.function = &get_expected_photons_root; + F.params = &pars; + gsl_root_fsolver_set(s, &F, Tmin, Tmax); do { @@ -5432,7 +5438,7 @@ int fit_event2(event *ev, double *xopt, double *fmin, int *id, size_t n, double nlopt_set_min_objective(opt,nlopt_nll2,&fpars); /* Guess the position and t0 of the event using the QUAD fitter. */ - status = quad(ev,pos,&t0,10000); + status = quad(ev,pos,&t0,10000,0.1); if (status || t0 > GTVALID || NORM(pos) > PSUP_RADIUS) { /* If the QUAD fitter fails or returns something outside the PSUP or @@ -5474,7 +5480,7 @@ int fit_event2(event *ev, double *xopt, double *fmin, int *id, size_t n, double ub[3] = GTVALID; /* Find the peaks in the Hough transform of the event. */ - find_peaks(ev, pos, 100, 100, peak_theta, peak_phi, &npeaks, LEN(peak_theta),0.1); + find_peaks(ev, pos, 100, 100, peak_theta, peak_phi, &npeaks, LEN(peak_theta), 0.1); /* Don't fit more than 3 peaks for now. */ if (npeaks > 3) npeaks = 3; |