From 73f6b9d059436497bb5f95d7b60919d18e120b6e Mon Sep 17 00:00:00 2001 From: tlatorre Date: Mon, 27 Aug 2018 11:20:26 -0500 Subject: update code to use get_index_snoman* functions to calculate the index of refraction --- src/likelihood.c | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) (limited to 'src/likelihood.c') diff --git a/src/likelihood.c b/src/likelihood.c index 4e24a79..a4864a8 100644 --- a/src/likelihood.c +++ b/src/likelihood.c @@ -78,7 +78,7 @@ static double gsl_muon_time(double x, void *params) { double *params2 = (double *) params; double T0 = params2[0]; - double pos0[3], dir[3], pos[3], pmt_dir[3]; + double pos0[3], dir[3], pos[3], pmt_dir[3], R, n, wavelength0; int i; double t; i = (int) params2[1]; @@ -93,13 +93,19 @@ static double gsl_muon_time(double x, void *params) pos[1] = pos0[1] + dir[1]*x; pos[2] = pos0[2] + dir[2]*x; + R = NORM(pos); + SUB(pmt_dir,pmts[i].pos,pos); double distance = NORM(pmt_dir); /* FIXME: I just calculate delta assuming 400 nm light. */ - double wavelength0 = 400.0; - double n = get_index(HEAVY_WATER_DENSITY, wavelength0, 10.0); + wavelength0 = 400.0; + if (R <= AV_RADIUS) { + n = get_index_snoman_d2o(wavelength0); + } else { + n = get_index_snoman_h2o(wavelength0); + } t = x/SPEED_OF_LIGHT + distance*n/SPEED_OF_LIGHT; @@ -197,7 +203,7 @@ double nll_muon(event *ev, double T, double *pos, double *dir, double t0) * just to split up the integral, we only need to find a point along * the track close enough such that the integral isn't completely zero. */ - theta_cerenkov = acos(1/get_index(WATER_DENSITY,400.0,10.0)); + theta_cerenkov = acos(1/get_index_snoman_d2o(400.0)); /* Now, we compute the distance along the track where the PMT is at the * Cerenkov angle. */ -- cgit