From b1076924108c7c9a109f803fcb989411e71f2695 Mon Sep 17 00:00:00 2001 From: Luke Pickering Date: Fri, 5 Apr 2019 10:17:08 -0500 Subject: [PATCH] Adds controls for energy range, regularization, and max off axis position. --- flux_hist.js | 78 ++++++++++++++++++++++++++++++++++++++++++++-------- osc.html | 60 +++++++++++++++++++++++++++++++++------- 2 files changed, 116 insertions(+), 22 deletions(-) diff --git a/flux_hist.js b/flux_hist.js index 75cb3fe..5efcdcc 100644 --- a/flux_hist.js +++ b/flux_hist.js @@ -72,6 +72,26 @@ class flux_hist { return points; } + + static FindBin(e_bins, E) { + console.log("FindBin: ", E); + if (E < e_bins[0]) { + return -1; + } + if (E >= e_bins[e_bins.length - 1]) { + return e_bins.length; + } + + for (let bi_it = 0; bi_it < e_bins.length; ++bi_it) { + if (E < e_bins[bi_it]) { + console.log("FindBin: bin = ", bi_it); + + return (bi_it - 1); + } + } + } + + FindBin(E) { return flux_hist.FindBin(this.e_bins, E); } }; class OffAxis_flux_hist { @@ -93,7 +113,7 @@ class OffAxis_flux_hist { this.mat = lalolib.transpose(lalolib.array2mat(bincontent)); } - flux_match(target) { + flux_match(target, regfac = 1E-9, emin = 0, emax = 5, oamax = 40) { if (target.bincontent.length != this.mat.m) { console.log(`Cannot perform flux match. ND EBin count = ${ @@ -101,29 +121,63 @@ class OffAxis_flux_hist { return; } + let lowbin = target.FindBin(emin); + let upbin = target.FindBin(emax); + let bigfac = 1E15; - let regfac = 1E-9; let targetvec = lalolib.array2vec(target.bincontent); let bigmat = lalolib.entrywisemul(this.mat, bigfac); targetvec = lalolib.entrywisemul(targetvec, bigfac); + + let emin_bin = 0; + let emax_bin = target.bincontent.length; + if (upbin < target.bincontent.length) { + emax_bin = upbin + 1; + } + if (lowbin > 0) { + emin_bin = lowbin; + } + + let upbin_oa = bigmat.n; + let oamax_bin = flux_hist.FindBin(this.oa_bins, oamax); + if (oamax_bin < (this.oa_bins.length - 1)) { + upbin_oa = oamax_bin + 1; + } + + console.log("ERange: ", emin, emax, "Bin range: ", lowbin, upbin, "max oa bin: ", oamax_bin, upbin_oa); + + console.log("Before slice: ", bigmat, targetvec); + + bigmat = lalolib.get(bigmat, lalolib.range(emin_bin, emax_bin), + lalolib.range(0, upbin_oa)); + targetvec = lalolib.get(targetvec, lalolib.range(emin_bin, emax_bin)); + + console.log("After slice: ", bigmat, targetvec); + let RHS = lalolib.mul(lalolib.transpose(bigmat), targetvec); // make the penalty matrix // there is probably a better way to do this - // e.g., off-diag + // e.g., off-diag let A = lalolib.eye(bigmat.n); - for ( let i = 0; i < A.m; i++ ) { - for ( let j = 0; j < A.n; j++ ) { - if ( i == j - 1 ) { - A[i, j] = -1; - } - } + for (let i = 0; i < A.m; i++) { + for (let j = 0; j < A.n; j++) { + if (i == j - 1) { + A[i, j] = -1; + } + } } - - A = lalolib.entrywisemul(A, bigfac*regfac); + + A = lalolib.entrywisemul(A, bigfac * regfac); let LHS = lalolib.add(lalolib.xtx(bigmat), lalolib.xtx(A)); let result = lalolib.solve(LHS, RHS); - return { bf: lalolib.mul(this.mat, result), coeffs: result }; + let smallmat = lalolib.entrywisemul(bigmat, 1.0 / bigfac); + + let bf = new flux_hist(target.e_bins.slice(emin_bin, emax_bin + 1), + lalolib.mul(smallmat, result)); + let coeffs = new flux_hist(this.oa_bins.slice(0, upbin_oa+1), + result); + return {bf : bf, coeffs : coeffs}; } }; diff --git a/osc.html b/osc.html index ffd4b72..b72eb9f 100644 --- a/osc.html +++ b/osc.html @@ -116,6 +116,10 @@ var pdg_from = 14; var pdg_to = 14; var baseline_km = 1300; + var regfac = 1E-9; + var emin_prism = 0; + var emax_prism = 5; + var maxoffaxis = 40; function AddCurve(plot, params, curve_state=0) { @@ -156,7 +160,7 @@ function OscProbInit(){ - op.InitializeTable(document.getElementById("param_table")); + // op.InitializeTable(document.getElementById("param_table")); let op_el = document.getElementById("OscProbPlot"); var op_plot = new OscProbPlot(); @@ -182,15 +186,16 @@ fp.Draw(FD_numode_numu_flux); FD_numode_numu_flux_osc.Oscillate(pdg_from, pdg_to, baseline_km, op); fp.Draw(FD_numode_numu_flux_osc); - let bfbinc = ND_numuode_numu_flux.flux_match(FD_numode_numu_flux_osc); - let bfflux = new flux_hist(DUNE_flux.FD.numode.numu.bins_e, bfbinc.bf); + let bfbinc = ND_numuode_numu_flux.flux_match(FD_numode_numu_flux_osc, regfac, emin_prism, emax_prism, maxoffaxis); + let bfflux = bfbinc.bf; bfflux.line_class = "orange_line"; + console.log(bfflux); fp.Draw(bfflux); let lincombplot_el = document.getElementById("LinCombPlot"); let lcp = new off_axis_lincomb_plot(); - let lincomb = new flux_hist(DUNE_flux.ND.numode.numu.bins_OA, bfbinc.coeffs); + let lincomb = bfbinc.coeffs; lincomb.line_class = "orange_line"; lcp.DrawAxes(lincombplot_el, lincomb, 1E5, 1); lcp.Draw(lincomb); @@ -202,17 +207,19 @@ FD_numode_numu_flux_osc.Oscillate(pdg_from, pdg_to, baseline_km, op); fp.Clear(2); fp.Draw(FD_numode_numu_flux_osc); - let bfbinc = ND_numuode_numu_flux.flux_match(FD_numode_numu_flux_osc); - let bfflux = new flux_hist(DUNE_flux.FD.numode.numu.bins_e, bfbinc.bf); + console.log("Min/Max: ", emin_prism, emax_prism); + let bfbinc = ND_numuode_numu_flux.flux_match(FD_numode_numu_flux_osc, regfac, emin_prism, emax_prism, maxoffaxis); + let bfflux = bfbinc.bf; bfflux.line_class = "orange_line"; fp.Draw(bfflux); lcp.Clear(1); - lincomb = new flux_hist(DUNE_flux.ND.numode.numu.bins_OA, bfbinc.coeffs); + let lincomb = bfbinc.coeffs; lincomb.line_class = "orange_line"; lcp.Draw(lincomb); }; + let cc = document.getElementById("ConstraintControls"); InitializeConstraintPlots(cc, function(chg_par) { op.Set(chg_par.xax_name, chg_par.xax_point); @@ -266,6 +273,23 @@ baseline_km = chosen; AddCurve(op_plot,op,1); }); + d3.select("#regfac_input").on('change',function(){ + let chosen = parseFloat(this.value); + regfac = chosen; + }); + d3.select("#emin_input").on('change',function(){ + let chosen = parseFloat(this.value); + emin_prism = chosen; + }); + d3.select("#emax_input").on('change',function(){ + let chosen = parseFloat(this.value); + emax_prism = chosen; + }); + d3.select("#maxoffaxis_input").on('change',function(){ + let chosen = parseFloat(this.value); + maxoffaxis = chosen; + }); + } @@ -326,7 +350,7 @@

Oscillation Probability

- +
@@ -350,12 +374,28 @@
Click/Drag to choose parameters
+
+ + +
+
+ + +
+
+ + +
+
+ + +
-
+