aboutsummaryrefslogtreecommitdiff
path: root/src/muon.c
AgeCommit message (Collapse)Author
2019-11-05update guess_energy()tlatorre
This commit updates guess_energy() which is used to seed the energy for the likelihood fit. Previously we estimated the energy by summing up the charge in a 42 degree cone around the proposed direction and then dividing that by 6 (since electrons in SNO and SNO+ produce approximately 6 hits/MeV). Now, guess_energy() estimates the energy by calculating the expected number of photons produced from Cerenkov light, EM showers, and delta rays for a given particle at a given energy. The most likely energy is found by bisecting the difference between the expected number of photons and the observed charge to find when they are equal. This improves things dramatically for cosmic muons which have energies of ~200 GeV. Previously the initial guess was always very low (~1 GeV) and the fit could take > 1 hour to increase the energy.
2019-09-24update shower position distribution parameters for muonstlatorre
This commit updates the a and b parameters for the gamma distribution used to describe the position distribution of shower photons produced along the direction of the muon. Previously I had been assuming b was equal to the radiation length and using a formula from the PDG to calculate a from that. However, this formula doesn't seem to be valid for muons (the formula comes from a section describing the shower profile of electrons and gammas, so it's not surprising). Therefore, now we don't assume any relationship between a and b. Now, the value of a is approximated by a constant since I couldn't find any functional relationship as a function of energy to describe a very well (and it's approximately constant), and b is approximated by a single degree polynomial fit to the values I got from simulating muons in RAT-PAC as a function of energy. Note that looking at the simulation data it seems like the position distribution of shower photons from muons isn't even very well described by a gamma distribution, so in the future it might be a good idea to come up with a better parameterization. Even if I stick with the gamma distribution, it would be good to revisit this in the future and fit for a and b over a wider range of energies.
2019-06-14set the maximum kinetic energy in the fit dynamically based on particle IDtlatorre
The range and energy loss tables have different maximum values for electrons, muons, and protons so we have to dynamically set the maximum energy of the fit in order to avoid a GSL interpolation error. This commit adds {electron,muon,proton}_get_max_energy() functions to return the maximum energy in the tables and that is then used to set the maximum value in the fit.
2019-05-22add a function to compute the number of shower photons for muonstlatorre
Similarly to electrons, I fit an analytic form to the ratio of the number of photons produced via shower particles over the radiative energy loss. In this case, I chose the functional form: ratio = a*(1-exp(-T/b)) since the ratio seemed to reach a constant value after a certain energy. I then simulated a 10 GeV muon and it appears that the ratio might actually decrease after that, so for higher energies I may have to come up with a different fit function.
2019-03-16add GPLv3 licensetlatorre
2019-03-07update code to allow you to run the fit outside of the src directorytlatorre
To enable the fitter to run outside of the src directory, I created a new function open_file() which works exactly like fopen() except that it searches for the file in both the current working directory and the path specified by an environment variable.
2019-01-31small updates to make sure we don't calculate nanstlatorre
2019-01-27add photons from delta rays to likelihood calculationtlatorre
This commit updates the likelihood function to take into account Cerenkov light produced from delta rays produced by muons. The angular distribution of this light is currently assumed to be constant along the track and parameterized in the same way as the Cerenkov light from an electromagnetic shower. Currently I assume the light is produced uniformly along the track which isn't exactly correct, but should be good enough.
2018-11-11update likelihood function to fit electrons!tlatorre
To characterize the angular distribution of photons from an electromagnetic shower I came up with the following functional form: f(cos_theta) ~ exp(-abs(cos_theta-mu)^alpha/beta) and fit this to data simulated using RAT-PAC at several different energies. I then fit the alpha and beta coefficients as a function of energy to the functional form: alpha = c0 + c1/log(c2*T0 + c3) beta = c0 + c1/log(c2*T0 + c3). where T0 is the initial energy of the electron in MeV and c0, c1, c2, and c3 are parameters which I fit. The longitudinal distribution of the photons generated from an electromagnetic shower is described by a gamma distribution: f(x) = x**(a-1)*exp(-x/b)/(Gamma(a)*b**a). This parameterization comes from the PDG "Passage of particles through matter" section 32.5. I also fit the data from my RAT-PAC simulation, but currently I am not using it, and instead using a simpler form to calculate the coefficients from the PDG (although I estimated the b parameter from the RAT-PAC data). I also sped up the calculation of the solid angle by making a lookup table since it was taking a significant fraction of the time to compute the likelihood function.
2018-10-18update fit to fit for electrons and protonstlatorre
2018-10-03add \ntlatorre
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-09-20delete unused variable distancetlatorre
2018-09-20add absorption lengths for D2O and H2O weighted by the Cerenkov spectrum and ↵tlatorre
the quantum efficiency
2018-09-18speed likelihood calculation up a bittlatorre
2018-09-17update likelihood to calculate absorption length correctlytlatorre
2018-09-17update muon kinetic energy calculationtlatorre
This commit updates the calculation of the muon kinetic energy as a function of distance along the track. Previously I was using an approximation from the PDG, but it doesn't seem to be very accurate and won't generalize to the case of electrons. The kinetic energy is now calculated using the tabulated values of dE/dx as a function of energy.
2018-09-17fix a bug in get_dEdx()tlatorre
2018-09-11add absorption lengthtlatorre
This commit adds the absorption length to the likelihood calculation. For now I'm just using a single number independent of wavelength. I should update this in the future to actually use the absorption lengths as measured by SNO and then calculate an overall absorption length weighted by the Cerenkov spectrum and the PMT quantum efficiency.
2018-09-11add PMT responsetlatorre
This commit adds code to read in the PMT response from the PMTR bank from SNOMAN. This file was used for the grey disk model in SNOMAN and was created using a full 3D simulation of the PMT and concentrator. Since the PMT response in SNOMAN included the quantum efficiency of the PMT, we have to divide that out to get just the PMT response independent of the quantum efficiency. I also updated the likelihood calculation to use the pmt response. Currently the energy is being fit too high which I think will improve when we update the solid angle calculation to use the radius of the concentrator instead of the PMT.
2018-08-31add muon critical energy for D2Otlatorre
2018-08-28add path to the likelihood fittlatorre
This commit updates the likelihood fit to use the KL path expansion. Currently, I'm just using one coefficient for the path in both x and y.
2018-08-27fix how multiple Coulomb scattering is treatedtlatorre
Previously I had been assuming that a particle undergoing many small angle Coulomb scatters had a track direction whose polar angle was a Gaussian. However, this was just due to a misunderstanding of the PDG section "Multiple scattering through small angles" in the "Passage of particles through matter" article. In fact, what is described by a Gaussian is the polar angle projected onto a plane. Therefore the distribution of the polar angle is actually: (1/(sqrt(2*pi)*theta0**2))*theta*exp(-theta**2/(2*theta0)) This commit updates the code in scattering.c to correctly calculate the probability that a photon is emitted at a particular angle. I also updated test-likelihood.c to simulate a track correctly.
2018-08-14fix how the RMS scattering angle is calculatedtlatorre
The RMS scattering angle calculation comes from Equation 33.15 in the PDG article on the passage of particles through matter. It's not entirely obvious if this equation is correct for a long track. It seems like it should be integrated along the track to add up the contributions at different energies, but it's not obvious how to do that with the log term. In any case, the way I was previously calculating it (by using the momentum and velocity at each point along the track) was definitely wrong. I will try this out and perhaps try to integrate it later.
2018-08-14add refractive index for heavy and light water from snomantlatorre
2018-08-14update muon.c to use the new table for heavy watertlatorre
2018-08-14move everything to src directorytlatorre