aboutsummaryrefslogtreecommitdiff
path: root/src
AgeCommit message (Collapse)Author
2018-10-21speed up get_total_charge_approx() by precomputing some variablestlatorre
2018-10-21fix use of uninitialized variablestlatorre
2018-10-19don't call path_init() when doing the fast likelihood calculation to speed ↵tlatorre
things up
2018-10-19add MIN_RATIO_FAST to speed up the "fast" likelihood calculationtlatorre
2018-10-19speed up get_total_charge_approx()tlatorre
2018-10-19epsrel -> npointstlatorre
2018-10-19add interp2d() for fast bilinear 2D interpolationtlatorre
2018-10-19speed up sin_theta calculation by replacing sin(acos(cos_theta)) with ↵tlatorre
sqrt(1-pow(cos_theta,2))
2018-10-19update path integral to use a fixed number of pointstlatorre
I noticed when fitting electrons that the cquad integration routine was not very stable, i.e. it would return different results for *very* small changes in the fit parameters which would cause the fit to stall. Since it's very important for the minimizer that the likelihood function not jump around, I am switching to integrating over the path by just using a fixed number of points and using the trapezoidal rule. This seems to be a lot more stable, and as a bonus I was able to combine the three integrals (direct charge, indirect charge, and time) so that we only have to do a single loop. This should hopefully make the speed comparable since the cquad routine was fairly effective at only using as many function evaluations as needed. Another benefit to this approach is that if needed, it will be easier to port to a GPU.
2018-10-18fix a bug in get_total_charge_approx()tlatorre
This commit fixes a bug which was double counting the pmt response when computing the direct charge and incorrectly multiplying the reflected charge by the pmt response. I think this was just a typo left in when I added the reflected charge.
2018-10-18update path_init() to check for a divide by zerotlatorre
This commit updates path_init() to check that beta > 0 before dividing by it to compute the time. Previously when fitting electrons it would occasionally divide by zero which would cause the inf to propagate all the way through the likelihood function.
2018-10-18make sure that the kinetic energy is zero at the last steptlatorre
Occasionally when fitting electrons the kinetic energy at the last step would be high enough that the electron never crossed the BETA_MIN threshold which would cause the gsl routine to throw an error. This commit updates particle_init() to set the kinetic energy at the last step to zero to make sure that we can bisect the point along the track where the speed drops to BETA_MIN.
2018-10-18hardcode the density when computing dE/dxtlatorre
Since we only have the range and dE/dx tables for light water for electrons and protons it's not correct to use the heavy water density. Also, even though we have both tables for muons, currently we only load the heavy water table, so we hardcode the density to that of heavy water. In the future, it would be nice to load both tables and use the correct one depending on if we are fitting in the heavy or light water.
2018-10-18fix the likelihood function to return the *negative* log likelihood of the ↵tlatorre
path coefficients Previously I was adding the log likelihood of the path coefficients instead of the *negative* log likelihood! When fitting electrons this would sometimes cause the fit to become unstable and continue increasing the path coefficients without bound since the gain in the likelihood caused by increasing the coefficients was more than the loss caused by a worse fit to the PMT data. Doh!
2018-10-18update theta0 calculation to use the radiation length in light watertlatorre
Previously I was using the radiation length in light water but scaling it by the density of heavy water, which isn't correct. Since the radiation length in heavy and light water is almost identical, we just use the radiation length in light water.
2018-10-18update fit to fit for electrons and protonstlatorre
2018-10-17fix a bug in the theta0 calculation for a pathtlatorre
This commit fixes a bug in the calculation of the average rms width of the angular distribution for a path with a KL expansion. I also made a lot of updates to the test-path program: - plot the distribution of the KL expansion coefficients - plot the standard deviation of the angular distribution as a function of distance along with the prediction - plot the simulated and reconstructed path in 3D
2018-10-12skip PMTs which weren't hit for the fast likelihood calculationtlatorre
2018-10-06add media_qoca_d2o_20060110.cmdtlatorre
2018-10-06prune output of mcvx and mctk bankstlatorre
2018-10-05epsrel = 1e-2 to speed up fittlatorre
2018-10-03update absorption coefficients from media_qoca_d2o_20060110.cmdtlatorre
2018-10-03add \ntlatorre
2018-10-03move python scripts into utils/ directorytlatorre
2018-10-03add a script to plot fit resultstlatorre
2018-10-03make sure wavelengths are uniformly spaced since we're using interp1d()tlatorre
2018-10-02update get_quantum_efficiency() to use interp1d()tlatorre
This speeds up the "boot up" time from ~30 seconds to ~12 seconds.
2018-10-02update MIN_THETA0 to 0.02tlatorre
2018-10-02update charge fraction to 0.5tlatorre
2018-10-02update optics to use LETA constants for the D2O phasetlatorre
This commit updates the PMT response bank and absorption lengths for H2O and D2O based on the values in mcprod which were used for the LETA analysis.
2018-10-01flag PMTs with the KPF_DIS bit settlatorre
Doh! I was previously using the tube number variable!
2018-10-01fall back to old scattering rms calculation when n = 0tlatorre
2018-10-01update negative log likelihood for path coefficientstlatorre
2018-10-01loop over all normal PMTs when calculating the expected number of photonstlatorre
Previously we ignored PMTs which were flagged when computing the expected number of PE for each PMT, but since we calculate the amount of reflected light here we need to include even PMTs which are offline (since they still reflect light).
2018-10-01use the PMT response table to calculate the amount of reflected lighttlatorre
To calculate the expected number of photons from reflected light we now integrate over the track and use the PMT response table to calculate what fraction of the light is reflected. Previously we were just using a constant fraction of the total detected light which was faster since we only had to integrate over the track once, but this should be more accurate.
2018-10-01add absorption length for acrylictlatorre
2018-10-01fix a bug when computing the absorption length in H2O and D2Otlatorre
Previously I was interpolating the absorption lengths using interp1d() but that only works when the x array is uniform. Since the wavelengths are not spaced uniformly, we have to use the GSL interpolation routines.
2018-09-26speed up fast likelihood calculationtlatorre
This commit updates the fast likelihood calculation to use the identity sin(a-b) = sin(a)*cos(b) - cos(a)*sin(b) to speed up the fast likelihood calculation.
2018-09-26speed up fast likelihood calculationtlatorre
This commit speeds up the fast likelihood calculation by avoiding calls to trigonometric functions where possible. Specifically we calculate sin(a) = sqrt(1-pow(cos(a),2)); instead of sin(a) = sin(acos(cos(a)));
2018-09-26increase number of points in cos(theta) interpolation to 1000tlatorre
2018-09-25update indirect scattering PDF start timetlatorre
Currently the PDF for scattered light is modelled as a flat distribution starting at some time t. Previously I was using the mean hit time for all PMTs, however this should really be a flat distribution in the time *residual* after the main peak. Therefore, the PDF now starts at the estimated time for direct photons.
2018-09-25update likelihood calculation to use PMT_TTS macrotlatorre
I accidentally hardcoded the single PE TTS to 1.5 ns in the likelihood calculation.
2018-09-25update integration bounds in likelihood calculationtlatorre
This commit updates the bounds of the track integration in the likelihood function to integrate up to 1 meter around the point at which the PMT is at the Cerenkov angle from the track. This fixes an issue I was seeing where a *very* small change in the fit paramters would cause the likelihood to jump by a large amount. I eventually tracked it down to the same issue I was seeing before which I solved by splitting up the integration into two intervals. However that fix did not seem to completely fix the issue. Based on initial tests with 500 MeV muons, this fix seems to do a much better job.
2018-09-25increase maxeval to 20 for the "quick" minimization phasetlatorre
When testing out the fitter on 500 MeV muons, there was at least one event which started out at a position very far from the true position. This event had a secondary electron like ring which is what I suspect caused the fit to start out in a position far from the true position. This fix correctly starts the minimization close to the true position. In the future I should look at updating get_direction() so that it finds the largest ring direction instead of just doing a weighted average of all the vectors from the position to the PMTs.
2018-09-21update likelihood function to include the probability of the path coefficientstlatorre
2018-09-21split up the track integral into two intervalstlatorre
This commit updates the likelihood calculation to split up the track integral into two intervals in some cases. I noticed when fitting some events that the likelihood value would change drastically for a very small change in the fit parameters. I eventually tracked it down to the fact that the track integral was occasionally returning a very small charge for a PMT which should have a very high charge. This was happening because the region of the track which was hitting the PMT was very small and the cquad integration routine was completely skipping it. The solution to this problem is a bit of a hack, but it seems to work. I first calculate where along the track (for a straight track) the PMT would be at the Cerenkov angle from the track. If this point is somewhere along the track then we split up the integral into two intervals: one going from the start of the track to this point and the other from the point to the end of the track. Since the cquad routine always samples points near the end of the intervals this should prevent it from completely skipping over the point in the track where the integrand is non-zero.
2018-09-21increase MAX_PE to 10000tlatorre
2018-09-20delete unused variable distancetlatorre
2018-09-20don't include the OWL PMTs in the likelihood calculationtlatorre
For some reason the OWL tubes have 9999.00 for the x, y, and z coordinates of the normal vector in the PMT file. For now, I'm just going to remove them from the likelihood calculation.
2018-09-20make sure direction vector is normalized in path_eval()tlatorre