From 15d7c926f81534edc107f2523ae2ff9fbe98a6e6 Mon Sep 17 00:00:00 2001 From: tlatorre Date: Mon, 7 Sep 2020 13:50:05 -0500 Subject: update a comment in chi2 --- utils/chi2 | 45 ++++++++++++++------------------------------- 1 file changed, 14 insertions(+), 31 deletions(-) diff --git a/utils/chi2 b/utils/chi2 index b3875d1..7a8342d 100755 --- a/utils/chi2 +++ b/utils/chi2 @@ -224,40 +224,23 @@ def make_nll(data, muons, mc, bins, print_nll=False): if any(x[i] < 0 for i in range(len(x))): return np.inf + # Get the data histograms. We need to do this within the + # likelihood function since we apply the energy bias corrections to the + # data. + # + # Note: Applying the energy bias correction to the data has a few + # issues. If the energy scale parameter changes by even a very small + # amount it can cause events to cross a bin boundary you will see a + # discrete jump in the likelihood. This isn't good for minimizers, and + # so when doing the first minimization with nlopt, I fix these + # parameters. Later, when we run the MCMC, these paramters are then + # free to float and the small discontinuities in the likelihood + # shouldn't be a problem. data_hists = get_data_hists_fast(ke_dict_data,x,bins) # Get the Monte Carlo histograms. We need to do this within the - # likelihood function since several of the parameters like the energy - # bias and resolution affect the histograms. - # - # Note: I really should be applying the bias term to the data instead - # of the Monte Carlo. The reason being that the reconstruction was all - # tuned based on the data from the Monte Carlo. For example, the single - # PE charge distribution is taken from the Monte Carlo. Therefore, if - # there is some difference in bias between data and Monte Carlo, we - # should apply it to the data. However, the reason I apply it to the - # Monte Carlo instead is because applying an energy bias correction to - # an analysis in which we're using histograms is not really a good - # idea. If the energy scale parameter changes by a very small amount - # and causes a single event to cross a bin boundary you will see a - # discrete jump in the likelihood. This isn't good for minimizers - # (although it may be ok with the MCMC depending on the size of the - # jump). To reduce this issue, I apply the energy bias correction to - # the Monte Carlo, and since there are almost 100 times more events in - # the Monte Carlo than in data, the issue is much smaller. - # - # All that being said, this doesn't completely get rid of the issue, - # but I do two things that should make it OK: - # - # 1. I use the SBPLX minimizer instead of COBYLA which should have a - # better chance of jumping over small discontinuities. - # - # 2. I ultimately run an MCMC and so if we get stuck somewhere close to - # the minimum, the MCMC results should correctly deal with any small - # discontinuities. - # - # Also, it's critical that I first adjust the data energy by whatever - # amount I find with the stopping muons and Michel distributions. + # likelihood function since we apply the energy resolution parameters + # to the Monte Carlo. mc_hists = get_mc_hists_fast(ke_dict,x,bins,apply_norm=True) muon_hists = get_mc_hists_fast(ke_dict_muon,x,bins,apply_norm=False) -- cgit