diff options
| -rwxr-xr-x | fr.py | 43 | ||||
| -rw-r--r-- | utils/fr.py | 8 | ||||
| -rw-r--r-- | utils/likelihood.py | 16 |
3 files changed, 37 insertions, 30 deletions
@@ -29,31 +29,32 @@ from utils.param import Param, ParamSet, get_paramsets def define_nuisance(): """Define the nuisance parameters.""" tag = ParamTag.SM_ANGLES + nuisance = [] g_prior = PriorsCateg.GAUSSIAN hg_prior = PriorsCateg.HALFGAUSS e = 1e-9 - nuisance = [ - Param(name='s_12_2', value=0.307, seed=[0.26, 0.35], ranges=[0., 1.], std=0.013, tex=r's_{12}^2', prior=g_prior, tag=tag), - Param(name='c_13_4', value=1-(0.02206)**2, seed=[0.995, 1-e], ranges=[0., 1.], std=0.00147, tex=r'c_{13}^4', prior=hg_prior, tag=tag), - Param(name='s_23_2', value=0.538, seed=[0.31, 0.75], ranges=[0., 1.], std=0.069, tex=r's_{23}^2', prior=g_prior, tag=tag), - Param(name='dcp', value=4.08404, seed=[0+e, 2*np.pi-e], ranges=[0., 2*np.pi], std=2.0, tex=r'\delta_{CP}', tag=tag), - Param( - name='m21_2', value=7.40E-23, seed=[7.2E-23, 7.6E-23], ranges=[6.80E-23, 8.02E-23], - std=2.1E-24, tex=r'\Delta m_{21}^2{\rm GeV}^{-2}', prior=g_prior, tag=tag - ), - Param( - name='m3x_2', value=2.494E-21, seed=[2.46E-21, 2.53E-21], ranges=[2.399E-21, 2.593E-21], - std=3.3E-23, tex=r'\Delta m_{3x}^2{\rm GeV}^{-2}', prior=g_prior, tag=tag - ) - ] + # nuisance.extend([ + # Param(name='s_12_2', value=0.307, seed=[0.26, 0.35], ranges=[0., 1.], std=0.013, tex=r's_{12}^2', prior=g_prior, tag=tag), + # Param(name='c_13_4', value=1-(0.02206)**2, seed=[0.995, 1-e], ranges=[0., 1.], std=0.00147, tex=r'c_{13}^4', prior=hg_prior, tag=tag), + # Param(name='s_23_2', value=0.538, seed=[0.31, 0.75], ranges=[0., 1.], std=0.069, tex=r's_{23}^2', prior=g_prior, tag=tag), + # Param(name='dcp', value=4.08404, seed=[0+e, 2*np.pi-e], ranges=[0., 2*np.pi], std=2.0, tex=r'\delta_{CP}', tag=tag), + # Param( + # name='m21_2', value=7.40E-23, seed=[7.2E-23, 7.6E-23], ranges=[6.80E-23, 8.02E-23], + # std=2.1E-24, tex=r'\Delta m_{21}^2{\rm GeV}^{-2}', prior=g_prior, tag=tag + # ), + # Param( + # name='m3x_2', value=2.494E-21, seed=[2.46E-21, 2.53E-21], ranges=[2.399E-21, 2.593E-21], + # std=3.3E-23, tex=r'\Delta m_{3x}^2{\rm GeV}^{-2}', prior=g_prior, tag=tag + # ) + # ]) tag = ParamTag.NUISANCE - nuisance.extend([ - Param(name='convNorm', value=1., seed=[0.5, 2. ], ranges=[0. , 50.], std=0.3, tag=tag), - Param(name='promptNorm', value=0., seed=[0. , 6. ], ranges=[0. , 50.], std=0.05, tag=tag), - Param(name='muonNorm', value=1., seed=[0.1, 2. ], ranges=[0. , 50.], std=0.1, tag=tag), - Param(name='astroNorm', value=6.9, seed=[0.1, 10.], ranges=[0. , 50.], std=0.1, tag=tag), - Param(name='astroDeltaGamma', value=2.5, seed=[1. , 3. ], ranges=[-5., 5. ], std=0.1, tag=tag) - ]) + # nuisance.extend([ + # Param(name='convNorm', value=1., seed=[0.5, 2. ], ranges=[0. , 50.], std=0.3, tag=tag), + # Param(name='promptNorm', value=0., seed=[0. , 6. ], ranges=[0. , 50.], std=0.05, tag=tag), + # Param(name='muonNorm', value=1., seed=[0.1, 2. ], ranges=[0. , 50.], std=0.1, tag=tag), + # Param(name='astroNorm', value=6.9, seed=[0.1, 10.], ranges=[0. , 50.], std=0.1, tag=tag), + # Param(name='astroDeltaGamma', value=2.5, seed=[1. , 3. ], ranges=[-5., 5. ], std=0.1, tag=tag) + # ]) return ParamSet(nuisance) diff --git a/utils/fr.py b/utils/fr.py index 51d5b6a..f529796 100644 --- a/utils/fr.py +++ b/utils/fr.py @@ -194,12 +194,13 @@ def normalise_fr(fr): array([ 0.16666667, 0.33333333, 0.5 ]) """ - return np.array(fr) / np.sum(fr) + return np.array(fr) / float(np.sum(fr)) def estimate_scale(args): """Estimate the scale at which new physics will enter.""" - m_eign = args.m3x_2 + try: m_eign = args.m3x_2 + except: m_eign = MASS_EIGENVALUES[1] if args.energy_dependance is EnergyDependance.MONO: scale = np.power( 10, np.round(np.log10(m_eign/args.energy)) - \ @@ -459,8 +460,7 @@ def u_to_fr(source_fr, matrix): """ composition = np.einsum( - 'ai, bi, a -> b', - np.abs(matrix)**2, np.abs(matrix)**2, source_fr, + 'ai, bi, a -> b', np.abs(matrix)**2, np.abs(matrix)**2, source_fr, ) ratio = composition / np.sum(source_fr) return ratio diff --git a/utils/likelihood.py b/utils/likelihood.py index cd1ead8..2e4a22d 100644 --- a/utils/likelihood.py +++ b/utils/likelihood.py @@ -115,12 +115,16 @@ def triangle_llh(theta, args, asimov_paramset, llh_paramset, fitter): ) m_eig_names = ['m21_2', 'm3x_2'] - mass_eigenvalues = [x.value for x in llh_paramset if x.name in m_eig_names] - ma_names = ['s_12_2', 'c_13_4', 's_23_2', 'dcp'] - sm_u = fr_utils.angles_to_u( - [x.value for x in llh_paramset if x.name in ma_names] - ) + + if set(m_eig_names+ma_names).issubset(set(llh_paramset.names)): + mass_eigenvalues = [x.value for x in llh_paramset if x.name in m_eig_names] + sm_u = fr_utils.angles_to_u( + [x.value for x in llh_paramset if x.name in ma_names] + ) + else: + mass_eigenvalues = fr_utils.MASS_EIGENVALUES + sm_u = fr_utils.NUFIT_U if args.energy_dependance is EnergyDependance.MONO: u = fr_utils.params_to_BSMu( @@ -164,6 +168,8 @@ def triangle_llh(theta, args, asimov_paramset, llh_paramset, fitter): param.value = flavour_angles[idx] print 'llh_paramset', llh_paramset + print 'hypo_paramset', hypo_paramset + print 'fr', fr if args.likelihood is Likelihood.FLAT: llh = 1. elif args.likelihood is Likelihood.GAUSSIAN: |