From 346ed95f507a6b99131f998867d49c20b6b0f6be Mon Sep 17 00:00:00 2001 From: Nicolas Tessore Date: Mon, 12 Jun 2023 10:35:14 +0100 Subject: [PATCH] BUG(core): fix extrapolation in trapz_product() (#105) Fixes the incorrect extrapolation in `glass.core.array.trapz_product()` by making sure the interval over which all functions are interpolated telescopes down to the union of the input intervals. This caused a bug in `glass.points.effective_bias()`. Fixes: #104 Fixed: Incorrect extrapolation in `glass.core.array.trapz_product()`, causing a bug in `glass.points.effective_bias()` --- glass/core/array.py | 3 ++- glass/core/test/test_array.py | 14 ++++++++++++++ 2 files changed, 16 insertions(+), 1 deletion(-) diff --git a/glass/core/array.py b/glass/core/array.py index 7d80f855..dd0e4cf5 100644 --- a/glass/core/array.py +++ b/glass/core/array.py @@ -53,7 +53,8 @@ def trapz_product(f, *ff, axis=-1): '''trapezoidal rule for a product of functions''' x, _ = f for x_, _ in ff: - x = np.union1d(x, x_[(x_ > x[0]) & (x_ < x[-1])]) + x = np.union1d(x[(x >= x_[0]) & (x <= x_[-1])], + x_[(x_ >= x[0]) & (x_ <= x[-1])]) y = np.interp(x, *f) for f_ in ff: y *= np.interp(x, *f_) diff --git a/glass/core/test/test_array.py b/glass/core/test/test_array.py index b5059bb7..348e0048 100644 --- a/glass/core/test/test_array.py +++ b/glass/core/test/test_array.py @@ -90,3 +90,17 @@ def test_ndinterp(): [[[2.15], [2.25], [2.35], [2.45]], [[2.45], [2.35], [2.25], [2.15]], [[2.15], [2.45], [2.25], [2.35]]]], atol=1e-15) + + +def test_trapz_product(): + from glass.core.array import trapz_product + + x1 = np.linspace(0, 2, 100) + f1 = np.full_like(x1, 2.0) + + x2 = np.linspace(1, 2, 10) + f2 = np.full_like(x2, 0.5) + + s = trapz_product((x1, f1), (x2, f2)) + + assert np.allclose(s, 1.0)