diff options
Diffstat (limited to 'src/likelihood.c')
| -rw-r--r-- | src/likelihood.c | 24 |
1 files changed, 24 insertions, 0 deletions
diff --git a/src/likelihood.c b/src/likelihood.c index 1bf9549..b06ebc4 100644 --- a/src/likelihood.c +++ b/src/likelihood.c @@ -70,6 +70,14 @@ particle *particle_init(int id, double T0, size_t n) dEdx = electron_get_dEdx(p->T[i-1], WATER_DENSITY); p->T[i] = p->T[i-1] - dEdx*(p->x[i]-p->x[i-1]); } + /* Make sure that the energy is zero at the last step. This is so that + * when we try to bisect the point along the track where the speed of + * the particle is equal to BETA_MIN, we guarantee that there is a + * point along the track where the speed drops below BETA_MIN. + * + * A possible future improvement would be to dynamically compute the + * range here using the dE/dx table instead of reading in the range. */ + p->T[n-1] = 0; break; case IDP_MU_MINUS: @@ -89,6 +97,14 @@ particle *particle_init(int id, double T0, size_t n) dEdx = muon_get_dEdx(p->T[i-1], HEAVY_WATER_DENSITY); p->T[i] = p->T[i-1] - dEdx*(p->x[i]-p->x[i-1]); } + /* Make sure that the energy is zero at the last step. This is so that + * when we try to bisect the point along the track where the speed of + * the particle is equal to BETA_MIN, we guarantee that there is a + * point along the track where the speed drops below BETA_MIN. + * + * A possible future improvement would be to dynamically compute the + * range here using the dE/dx table instead of reading in the range. */ + p->T[n-1] = 0; break; case IDP_PROTON: @@ -104,6 +120,14 @@ particle *particle_init(int id, double T0, size_t n) dEdx = proton_get_dEdx(p->T[i-1], WATER_DENSITY); p->T[i] = p->T[i-1] - dEdx*(p->x[i]-p->x[i-1]); } + /* Make sure that the energy is zero at the last step. This is so that + * when we try to bisect the point along the track where the speed of + * the particle is equal to BETA_MIN, we guarantee that there is a + * point along the track where the speed drops below BETA_MIN. + * + * A possible future improvement would be to dynamically compute the + * range here using the dE/dx table instead of reading in the range. */ + p->T[n-1] = 0; break; default: |