From 900a278192d09dee5e76ef7993083664cd551de2 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 12:02:55 +0200 Subject: [PATCH 01/10] Make keyword non-optional, fix docstrings for boost models --- clmm/utils/boost.py | 30 +++++++++++++++--------------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index 0abaa8a1f..174c45c6b 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -2,8 +2,8 @@ import numpy as np -def compute_nfw_boost(rvals, rscale=1000, boost0=0.1): - """Given a list of `rvals`, and optional `rscale` and `boost0`, return the corresponding +def compute_nfw_boost(rvals, rscale, boost_norm): + """Given a list of `rvals`, and optional `rscale` and `boost0`, returns the corresponding boost factor at each rval Parameters @@ -11,14 +11,14 @@ def compute_nfw_boost(rvals, rscale=1000, boost0=0.1): rvals : array_like Radii rscale : float, optional - scale radius for NFW in same units as rvals (default 2000 kpc) - boost0 : float, optional - Boost factor at each value of rvals + scale radius in same units as rvals + boost_norm : float, optional + Boost factor normalisation Returns ------- array - Boost factor + Boost factor at each value of rvals """ r_norm = np.array(rvals) / rscale @@ -39,7 +39,7 @@ def _calc_finternal(r_norm): return 1.0 + boost0 * (1 - _calc_finternal(r_norm)) / (r_norm**2 - 1) -def compute_powerlaw_boost(rvals, rscale=1000, boost0=0.1, alpha=-1.0): +def compute_powerlaw_boost(rvals, rscale, boost_norm, slope): """Given a list of `rvals`, and optional `rscale` and `boost0`, and `alpha`, return the corresponding boost factor at each `rval` @@ -47,22 +47,22 @@ def compute_powerlaw_boost(rvals, rscale=1000, boost0=0.1, alpha=-1.0): ---------- rvals : array_like Radii - rscale : float, optional - Scale radius for NFW in same units as rvals (default 2000 kpc) - boost0 : float, optional - Boost factor at each value of rvals - alpha : float, optional - Exponent from Melchior+16. Default: -1.0 + rscale : float + Scale radius in same units as rvals + boost_norm : float + Boost factor normalisation + slope : float + Exponent for the power-law parametrisation. NB: Melchior+16 uses -1.0 Returns ------- array - Boost factor + Boost factor at each value of rvals """ r_norm = np.array(rvals) / rscale - return 1.0 + boost0 * (r_norm) ** alpha + return 1.0 + boost0 * (r_norm) ** slope boost_models = { From b627ab752892387f5727382cb180198f94bdd4e6 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 12:50:09 +0200 Subject: [PATCH 02/10] Update to boost.py --- clmm/utils/boost.py | 46 ++++++++++++++++++++++++++++++--------------- 1 file changed, 31 insertions(+), 15 deletions(-) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index 174c45c6b..c427eba9e 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -2,17 +2,31 @@ import numpy as np -def compute_nfw_boost(rvals, rscale, boost_norm): - """Given a list of `rvals`, and optional `rscale` and `boost0`, returns the corresponding - boost factor at each rval +def compute_nfw_boost(rvals, rscale, boost0=0.1): + r"""Computes the boost factor at radii rvals using a parametric form + following that of the analytical NFW excess surface density. + + Setting :math:`x = R/R_s` + + .. math:: + B(x) = 1 + B_0\frac{1-F(x)}{(x)^2 -1} + + where + + .. math:: + F(x) = \begin{cases} \frac{\arctan \sqrt{x^2 -1 }}{\sqrt{x^2 -1 }} & \text{if $x>1$}, \\ + \text{1} & \text{if x = 1}, \\ + \frac{\text{arctanh} \sqrt{1-x^2 }}{\sqrt{1- x^2 }} & \text{if $x<1$}. + \end{cases} + Parameters ---------- rvals : array_like Radii - rscale : float, optional - scale radius in same units as rvals - boost_norm : float, optional + rscale : float + Scale radius in same units as rvals + boost0: float, optional Boost factor normalisation Returns @@ -39,9 +53,11 @@ def _calc_finternal(r_norm): return 1.0 + boost0 * (1 - _calc_finternal(r_norm)) / (r_norm**2 - 1) -def compute_powerlaw_boost(rvals, rscale, boost_norm, slope): - """Given a list of `rvals`, and optional `rscale` and `boost0`, and `alpha`, - return the corresponding boost factor at each `rval` +def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): + r"""Computes the boost factor at radii rvals using a power-law parametric form + + .. math:: + B(R) = 1 + B_0 \left(\frac{R}{R_s}\right)^\alpha Parameters ---------- @@ -49,10 +65,10 @@ def compute_powerlaw_boost(rvals, rscale, boost_norm, slope): Radii rscale : float Scale radius in same units as rvals - boost_norm : float + boost0 : float, optional Boost factor normalisation - slope : float - Exponent for the power-law parametrisation. NB: Melchior+16 uses -1.0 + alpha : float, optional + Exponent for the power-law parametrisation. NB: -1.0 (as in Melchior+16) Returns ------- @@ -62,7 +78,7 @@ def compute_powerlaw_boost(rvals, rscale, boost_norm, slope): r_norm = np.array(rvals) / rscale - return 1.0 + boost0 * (r_norm) ** slope + return 1.0 + boost0 * (r_norm) ** alpha boost_models = { @@ -91,7 +107,7 @@ def correct_sigma_with_boost_values(sigma_vals, boost_factors): return sigma_corrected -def correct_sigma_with_boost_model(rvals, sigma_vals, boost_model="nfw_boost", **boost_model_kw): +def correct_sigma_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale, **boost_model_kw): """Given a boost model and sigma profile, compute corrected sigma Parameters @@ -112,7 +128,7 @@ def correct_sigma_with_boost_model(rvals, sigma_vals, boost_model="nfw_boost", * correted radial profile """ boost_model_func = boost_models[boost_model] - boost_factors = boost_model_func(rvals, **boost_model_kw) + boost_factors = boost_model_func(rvals, boost_rscale, **boost_model_kw) sigma_corrected = np.array(sigma_vals) / boost_factors return sigma_corrected From e96248826e661822c571898b65b7285a06c5a052 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 12:50:36 +0200 Subject: [PATCH 03/10] Adapt notebook to new arguments/keywords --- examples/demo_boost_factors.ipynb | 59 ++++++++++++++++++++++++++----- 1 file changed, 51 insertions(+), 8 deletions(-) diff --git a/examples/demo_boost_factors.ipynb b/examples/demo_boost_factors.ipynb index 6f236bd22..5d5e460a9 100644 --- a/examples/demo_boost_factors.ipynb +++ b/examples/demo_boost_factors.ipynb @@ -2,6 +2,7 @@ "cells": [ { "cell_type": "markdown", + "id": "6afbf4b7", "metadata": {}, "source": [ "# Boost factors" @@ -10,6 +11,7 @@ { "cell_type": "code", "execution_count": null, + "id": "e3cd679b", "metadata": {}, "outputs": [], "source": [ @@ -24,6 +26,7 @@ }, { "cell_type": "markdown", + "id": "257a1892", "metadata": {}, "source": [ "Make sure we know which version we're using" @@ -32,6 +35,7 @@ { "cell_type": "code", "execution_count": null, + "id": "2ebc9e45", "metadata": {}, "outputs": [], "source": [ @@ -40,6 +44,7 @@ }, { "cell_type": "markdown", + "id": "67544ec1", "metadata": {}, "source": [ "### Define cosmology object" @@ -48,6 +53,7 @@ { "cell_type": "code", "execution_count": null, + "id": "cb23750b", "metadata": {}, "outputs": [], "source": [ @@ -56,6 +62,7 @@ }, { "cell_type": "markdown", + "id": "7b50a10e", "metadata": {}, "source": [ "First, we want to generate a $\\Delta\\Sigma$ (excess surface density) profile from mock data, to which we can apply boost factors. The mock data is generated in the following cells." @@ -63,6 +70,7 @@ }, { "cell_type": "markdown", + "id": "f4b319de", "metadata": {}, "source": [ "Generate cluster object from mock data" @@ -71,6 +79,7 @@ { "cell_type": "code", "execution_count": null, + "id": "2055117e", "metadata": {}, "outputs": [], "source": [ @@ -98,6 +107,7 @@ }, { "cell_type": "markdown", + "id": "8e24cea1", "metadata": {}, "source": [ "Loading this into a CLMM cluster object centered on (0,0)" @@ -106,6 +116,7 @@ { "cell_type": "code", "execution_count": null, + "id": "a62553b4", "metadata": {}, "outputs": [], "source": [ @@ -116,6 +127,7 @@ }, { "cell_type": "markdown", + "id": "a9f5c58f", "metadata": {}, "source": [ "Compute cross and tangential excess surface density for each source galaxy" @@ -124,6 +136,7 @@ { "cell_type": "code", "execution_count": null, + "id": "c961c355", "metadata": {}, "outputs": [], "source": [ @@ -141,6 +154,7 @@ }, { "cell_type": "markdown", + "id": "2aef9237", "metadata": {}, "source": [ "Calculate the binned profile" @@ -149,6 +163,7 @@ { "cell_type": "code", "execution_count": null, + "id": "8dbb5041", "metadata": {}, "outputs": [], "source": [ @@ -175,6 +190,7 @@ }, { "cell_type": "markdown", + "id": "d1e5e1d8", "metadata": {}, "source": [ "Plot the $\\Delta\\Sigma$ profile" @@ -183,6 +199,7 @@ { "cell_type": "code", "execution_count": null, + "id": "afb483c9", "metadata": {}, "outputs": [], "source": [ @@ -202,6 +219,7 @@ }, { "cell_type": "markdown", + "id": "8323cb3b", "metadata": {}, "source": [ "## Boost Factors" @@ -209,6 +227,7 @@ }, { "cell_type": "markdown", + "id": "de848b74", "metadata": {}, "source": [ "CLMM offers two boost models, the NFW boost model, and a powerlaw boost model. \n", @@ -223,6 +242,7 @@ { "cell_type": "code", "execution_count": null, + "id": "fab7bd2a", "metadata": {}, "outputs": [], "source": [ @@ -233,8 +253,19 @@ ")" ] }, + { + "cell_type": "code", + "execution_count": null, + "id": "d10a993d-acbc-4525-807d-9f2851bf48d8", + "metadata": {}, + "outputs": [], + "source": [ + "print(u.compute_nfw_boost(1000, 1000, boost0=0.1))" + ] + }, { "cell_type": "markdown", + "id": "2e356f74", "metadata": {}, "source": [ "Plot the two boost factors, $B(R)$" @@ -243,19 +274,22 @@ { "cell_type": "code", "execution_count": null, + "id": "05fd04e0", "metadata": {}, "outputs": [], "source": [ - "plt.plot(cl.DeltaSigma_profile[\"radius\"], nfw_boost, label=\"NFW boost factor\")\n", - "plt.plot(cl.DeltaSigma_profile[\"radius\"], powerlaw_boost, label=\"Powerlaw boost factor\")\n", + "plt.loglog(cl.DeltaSigma_profile[\"radius\"], nfw_boost, label=\"NFW boost factor\", marker='.')\n", + "plt.loglog(cl.DeltaSigma_profile[\"radius\"], powerlaw_boost, label=\"Powerlaw boost factor\")\n", "plt.xlabel(\"Radius [kpc]\")\n", "plt.ylabel(\"$B(R)$\")\n", + "plt.axvline(1000)\n", "plt.legend()\n", "plt.show()" ] }, { "cell_type": "markdown", + "id": "e48db29d", "metadata": {}, "source": [ "The $\\Delta\\Sigma$ profiles can be corrected with the boost factor using `correct_sigma_with_boost_values` or `correct_sigma_with_boost_model`. \n", @@ -266,6 +300,7 @@ }, { "cell_type": "markdown", + "id": "f6e65d7b", "metadata": {}, "source": [ "\n", @@ -285,6 +320,7 @@ }, { "cell_type": "markdown", + "id": "e57a94d2", "metadata": {}, "source": [ "First we will apply the boost factor with `correct_sigma_with_boost_values`" @@ -293,6 +329,7 @@ { "cell_type": "code", "execution_count": null, + "id": "90b9c8d2", "metadata": {}, "outputs": [], "source": [ @@ -306,6 +343,7 @@ }, { "cell_type": "markdown", + "id": "28040f51", "metadata": {}, "source": [ "Plot the result" @@ -314,6 +352,7 @@ { "cell_type": "code", "execution_count": null, + "id": "65fcd8b7", "metadata": {}, "outputs": [], "source": [ @@ -342,7 +381,7 @@ " color=\"grey\",\n", ")\n", "\n", - "# plt.loglog()\n", + "#plt.loglog()\n", "plt.title(\"DeltaSigma profile\")\n", "plt.xlabel(\"Radius [kpc]\")\n", "plt.ylabel(\"$\\Delta\\Sigma [M_\\odot\\; Mpc^{-2}]$\")\n", @@ -352,6 +391,7 @@ }, { "cell_type": "markdown", + "id": "67c5e14a", "metadata": {}, "source": [ "Now the same again but with `correct_sigma_with_boost_model`" @@ -360,18 +400,21 @@ { "cell_type": "code", "execution_count": null, + "id": "115a8687", "metadata": {}, "outputs": [], "source": [ "Sigma_corrected_powerlaw_boost = u.correct_sigma_with_boost_model(\n", " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", - " boost_model=\"powerlaw_boost\",\n", + " \"powerlaw_boost\", # boost_model\n", + " 1000, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", ")\n", "Sigma_corrected_nfw_boost = u.correct_sigma_with_boost_model(\n", " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", - " boost_model=\"nfw_boost\",\n", + " \"nfw_boost\", # boost_model\n", + " 1000, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", ")\n", "\n", "plt.errorbar(\n", @@ -410,9 +453,9 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "wrk", "language": "python", - "name": "python3" + "name": "wrk" }, "language_info": { "codemirror_mode": { @@ -424,7 +467,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.6.9" + "version": "3.10.9" } }, "nbformat": 4, From 71d1c5085f9af0cfbc5aba56fad13e30e06ee6bf Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 12:57:47 +0200 Subject: [PATCH 04/10] update test_utils accordingly --- tests/test_utils.py | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/tests/test_utils.py b/tests/test_utils.py index a70524472..a605b984b 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -23,8 +23,9 @@ def test_compute_nfw_boost(): """Test the nfw model for boost factor""" # Test data rvals = np.arange(1, 11).tolist() + rscale = 1000 - boost_factors = utils.compute_nfw_boost(rvals) + boost_factors = utils.compute_nfw_boost(rvals, rscale) test_boost_factors = np.array( [ @@ -49,8 +50,9 @@ def test_compute_powerlaw_boost(): """Test the powerlaw model for boost factor""" # Test data rvals = np.arange(1, 11).tolist() # Cannot contain 0 due to reciprocal term + rscale = 1000 - boost_factors = utils.compute_powerlaw_boost(rvals) + boost_factors = utils.compute_powerlaw_boost(rvals, rscale) test_boost_factors = np.array( [101.0, 51.0, 34.33333333, 26.0, 21.0, 17.66666667, 15.28571429, 13.5, 12.11111111, 11.0] @@ -77,17 +79,18 @@ def test_correct_sigma_with_boost_model(): # Make test data rvals = np.arange(1, 11).tolist() sigma_vals = (2 ** np.arange(10)).tolist() + boost_rscale = 1000 for boost_model in utils.boost_models.keys(): # Check for no nans or inf with positive-definite rvals and sigma vals assert np.all( np.isfinite( - utils.correct_sigma_with_boost_model(rvals, sigma_vals, boost_model=boost_model) + utils.correct_sigma_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale) ) ) # Test requesting unsupported boost model - assert_raises(KeyError, utils.correct_sigma_with_boost_model, rvals, sigma_vals, "glue") + assert_raises(KeyError, utils.correct_sigma_with_boost_model, rvals, sigma_vals, "glue", boost_rscale) def test_compute_radial_averages(): From bba860e04a8d352fec405afa0e71d354759eaba7 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 15:38:52 +0200 Subject: [PATCH 05/10] update boost correction function names --- clmm/utils/__init__.py | 4 ++-- clmm/utils/boost.py | 39 ++++++++++++++++--------------- examples/demo_boost_factors.ipynb | 10 ++++---- tests/test_utils.py | 12 ++++++---- 4 files changed, 34 insertions(+), 31 deletions(-) diff --git a/clmm/utils/__init__.py b/clmm/utils/__init__.py index bcd1806ac..a2c0854bb 100644 --- a/clmm/utils/__init__.py +++ b/clmm/utils/__init__.py @@ -13,8 +13,8 @@ from .boost import ( compute_nfw_boost, compute_powerlaw_boost, - correct_sigma_with_boost_values, - correct_sigma_with_boost_model, + correct_with_boost_values, + correct_with_boost_model, boost_models, ) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index c427eba9e..5fc702dec 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -1,4 +1,5 @@ """General utility functions that are used in multiple modules""" + import numpy as np @@ -54,8 +55,8 @@ def _calc_finternal(r_norm): def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): - r"""Computes the boost factor at radii rvals using a power-law parametric form - + r"""Computes the boost factor at radii rvals using a power-law parametric form + .. math:: B(R) = 1 + B_0 \left(\frac{R}{R_s}\right)^\alpha @@ -87,35 +88,35 @@ def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): } -def correct_sigma_with_boost_values(sigma_vals, boost_factors): - """Given a list of boost values and sigma profile, compute corrected sigma +def correct_with_boost_values(profile_vals, boost_factors): + """Given a list of profile (shear or DeltaSigma) values and boost values, compute corrected profile Parameters ---------- - sigma_vals : array_like - uncorrected sigma with cluster member dilution + profile_vals : array_like + Uncorrected profile values boost_factors : array_like - Boost values pre-computed + Boost values, pre-computed Returns ------- - sigma_corrected : numpy.ndarray - correted radial profile + profile_vals_corrected : numpy.ndarray + Corrected radial profile """ - sigma_corrected = np.array(sigma_vals) / np.array(boost_factors) - return sigma_corrected + profile_vals_corrected = np.array(profile_vals) / np.array(boost_factors) + return profile_vals_corrected -def correct_sigma_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale, **boost_model_kw): +def correct_with_boost_model(rvals, profile_vals, boost_model, boost_rscale, **boost_model_kw): """Given a boost model and sigma profile, compute corrected sigma Parameters ---------- rvals : array_like - radii - sigma_vals : array_like - uncorrected sigma with cluster member dilution + Radii + profile_vals : array_like + Uncorrected profile values boost_model : str, optional Boost model to use for correcting sigma @@ -124,14 +125,14 @@ def correct_sigma_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale, Returns ------- - sigma_corrected : numpy.ndarray - correted radial profile + profile_vals_corrected : numpy.ndarray + Correted radial profile """ boost_model_func = boost_models[boost_model] boost_factors = boost_model_func(rvals, boost_rscale, **boost_model_kw) - sigma_corrected = np.array(sigma_vals) / boost_factors - return sigma_corrected + profile_vals_corrected = np.array(profile_vals) / boost_factors + return profile_vals_corrected boost_models = { diff --git a/examples/demo_boost_factors.ipynb b/examples/demo_boost_factors.ipynb index 5d5e460a9..a7314ca7e 100644 --- a/examples/demo_boost_factors.ipynb +++ b/examples/demo_boost_factors.ipynb @@ -333,10 +333,10 @@ "metadata": {}, "outputs": [], "source": [ - "Sigma_corrected_powerlaw_boost = u.correct_sigma_with_boost_values(\n", + "Sigma_corrected_powerlaw_boost = u.correct_with_boost_values(\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"], powerlaw_boost\n", ")\n", - "Sigma_corrected_nfw_boost = u.correct_sigma_with_boost_values(\n", + "Sigma_corrected_nfw_boost = u.correct_with_boost_values(\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"], nfw_boost\n", ")" ] @@ -394,7 +394,7 @@ "id": "67c5e14a", "metadata": {}, "source": [ - "Now the same again but with `correct_sigma_with_boost_model`" + "Now the same again but with `correct_with_boost_model`" ] }, { @@ -404,13 +404,13 @@ "metadata": {}, "outputs": [], "source": [ - "Sigma_corrected_powerlaw_boost = u.correct_sigma_with_boost_model(\n", + "Sigma_corrected_powerlaw_boost = u.correct_with_boost_model(\n", " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", " \"powerlaw_boost\", # boost_model\n", " 1000, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", ")\n", - "Sigma_corrected_nfw_boost = u.correct_sigma_with_boost_model(\n", + "Sigma_corrected_nfw_boost = u.correct_with_boost_model(\n", " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", " \"nfw_boost\", # boost_model\n", diff --git a/tests/test_utils.py b/tests/test_utils.py index a605b984b..f4d41f2cb 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -62,7 +62,7 @@ def test_compute_powerlaw_boost(): assert_allclose(boost_factors, test_boost_factors) -def test_correct_sigma_with_boost_values(): +def test_correct_with_boost_values(): """ """ # Make test data rvals = np.arange(1, 11) @@ -70,11 +70,11 @@ def test_correct_sigma_with_boost_values(): test_unit_boost_factors = np.ones(rvals.shape).tolist() - corrected_sigma = utils.correct_sigma_with_boost_values(sigma_vals, test_unit_boost_factors) + corrected_sigma = utils.correct_with_boost_values(sigma_vals, test_unit_boost_factors) assert_allclose(sigma_vals, corrected_sigma) -def test_correct_sigma_with_boost_model(): +def test_correct_with_boost_model(): """ """ # Make test data rvals = np.arange(1, 11).tolist() @@ -85,12 +85,14 @@ def test_correct_sigma_with_boost_model(): # Check for no nans or inf with positive-definite rvals and sigma vals assert np.all( np.isfinite( - utils.correct_sigma_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale) + utils.correct_with_boost_model(rvals, sigma_vals, boost_model, boost_rscale) ) ) # Test requesting unsupported boost model - assert_raises(KeyError, utils.correct_sigma_with_boost_model, rvals, sigma_vals, "glue", boost_rscale) + assert_raises( + KeyError, utils.correct_with_boost_model, rvals, sigma_vals, "glue", boost_rscale + ) def test_compute_radial_averages(): From 5a894d4755a14ac2031c95ea728520ef1005240a Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 3 Oct 2024 15:50:35 +0200 Subject: [PATCH 06/10] pylint --- clmm/utils/boost.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index 5fc702dec..6a307981b 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -89,7 +89,8 @@ def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): def correct_with_boost_values(profile_vals, boost_factors): - """Given a list of profile (shear or DeltaSigma) values and boost values, compute corrected profile + """Given a list of profile values (e.g., shear or DeltaSigma) and boost values, + computes corrected profile Parameters ---------- From b91fc1d555e7005c5a7f44dd6abb78092e730a85 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 10 Oct 2024 15:38:09 +0200 Subject: [PATCH 07/10] update variable names and corresponding docstrings --- clmm/utils/boost.py | 44 +++++++++++++++++++++++--------------------- 1 file changed, 23 insertions(+), 21 deletions(-) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index 6a307981b..781fb287c 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -3,11 +3,11 @@ import numpy as np -def compute_nfw_boost(rvals, rscale, boost0=0.1): - r"""Computes the boost factor at radii rvals using a parametric form +def compute_nfw_boost(r_proj, r_scale, boost0=0.1): + r"""Computes the boost factor at radii r_proj using a parametric form following that of the analytical NFW excess surface density. - Setting :math:`x = R/R_s` + Setting :math:`x = r_{\rm proj}/r_{\rm scale}` .. math:: B(x) = 1 + B_0\frac{1-F(x)}{(x)^2 -1} @@ -23,20 +23,20 @@ def compute_nfw_boost(rvals, rscale, boost0=0.1): Parameters ---------- - rvals : array_like + r_proj : array_like Radii - rscale : float - Scale radius in same units as rvals + r_scale : float + Scale radius in same units as r_proj boost0: float, optional Boost factor normalisation Returns ------- array - Boost factor at each value of rvals + Boost factor at each value of r_proj """ - r_norm = np.array(rvals) / rscale + r_norm = np.array(r_proj) / r_scale def _calc_finternal(r_norm): radicand = r_norm**2 - 1 @@ -54,18 +54,18 @@ def _calc_finternal(r_norm): return 1.0 + boost0 * (1 - _calc_finternal(r_norm)) / (r_norm**2 - 1) -def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): - r"""Computes the boost factor at radii rvals using a power-law parametric form +def compute_powerlaw_boost(r_proj, r_scale, boost0=0.1, alpha=-1): + r"""Computes the boost factor at radii r_proj using a power-law parametric form .. math:: - B(R) = 1 + B_0 \left(\frac{R}{R_s}\right)^\alpha + B(r_{\rm proj}) = 1 + B_0 \left(\frac{r_{\rm proj}}{r_{\rm scale}}\right)^\alpha Parameters ---------- - rvals : array_like + r_proj : array_like Radii - rscale : float - Scale radius in same units as rvals + r_scale : float + Scale radius in same units as r_proj boost0 : float, optional Boost factor normalisation alpha : float, optional @@ -74,10 +74,10 @@ def compute_powerlaw_boost(rvals, rscale, boost0=0.1, alpha=-1): Returns ------- array - Boost factor at each value of rvals + Boost factor at each value of r_proj """ - r_norm = np.array(rvals) / rscale + r_norm = np.array(r_proj) / r_scale return 1.0 + boost0 * (r_norm) ** alpha @@ -109,20 +109,22 @@ def correct_with_boost_values(profile_vals, boost_factors): return profile_vals_corrected -def correct_with_boost_model(rvals, profile_vals, boost_model, boost_rscale, **boost_model_kw): +def correct_with_boost_model(r_proj, profile_vals, boost_model, boost_rscale, **boost_model_kw): """Given a boost model and sigma profile, compute corrected sigma Parameters ---------- - rvals : array_like + r_proj : array_like Radii profile_vals : array_like Uncorrected profile values - boost_model : str, optional + boost_model : str Boost model to use for correcting sigma - * 'nfw_boost' - NFW profile model (Default) * 'powerlaw_boost' - Powerlaw profile + boost_rscale : float + Scale radius to be used with the boost model + Returns ------- @@ -130,7 +132,7 @@ def correct_with_boost_model(rvals, profile_vals, boost_model, boost_rscale, **b Correted radial profile """ boost_model_func = boost_models[boost_model] - boost_factors = boost_model_func(rvals, boost_rscale, **boost_model_kw) + boost_factors = boost_model_func(r_proj, boost_rscale, **boost_model_kw) profile_vals_corrected = np.array(profile_vals) / boost_factors return profile_vals_corrected From cce412eb05a8649d69736e64c858d4bab67e80a7 Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 10 Oct 2024 15:46:06 +0200 Subject: [PATCH 08/10] fix the case r_proj=r_s in boost computation --- clmm/utils/boost.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/clmm/utils/boost.py b/clmm/utils/boost.py index 781fb287c..78e22c695 100644 --- a/clmm/utils/boost.py +++ b/clmm/utils/boost.py @@ -49,9 +49,11 @@ def _calc_finternal(r_norm): / (2 * np.lib.scimath.sqrt(radicand)) ) - return np.nan_to_num(finternal, copy=False, nan=1.0).real + return np.nan_to_num((1 - finternal) / radicand, copy=False, nan=1.0 / 3.0).real + # return np.nan_to_num(finternal, copy=False, nan=1.0).real - return 1.0 + boost0 * (1 - _calc_finternal(r_norm)) / (r_norm**2 - 1) + return 1.0 + boost0 * _calc_finternal(r_norm) + # return 1.0 + boost0 * (1 - _calc_finternal(r_norm)) / (r_norm**2 - 1) def compute_powerlaw_boost(r_proj, r_scale, boost0=0.1, alpha=-1): From 57d96d077093eee4f6b33c210ad30f42dc5d6b4d Mon Sep 17 00:00:00 2001 From: Celine Combet Date: Thu, 10 Oct 2024 15:51:28 +0200 Subject: [PATCH 09/10] Small update to notebook --- examples/demo_boost_factors.ipynb | 29 ++++++++++------------------- 1 file changed, 10 insertions(+), 19 deletions(-) diff --git a/examples/demo_boost_factors.ipynb b/examples/demo_boost_factors.ipynb index a7314ca7e..d7f8b320e 100644 --- a/examples/demo_boost_factors.ipynb +++ b/examples/demo_boost_factors.ipynb @@ -232,7 +232,8 @@ "source": [ "CLMM offers two boost models, the NFW boost model, and a powerlaw boost model. \n", "\n", - "Note that `compute_nfw_boost` requires two parameters to be specified, `rs` and `b0`, and `compute_powerlaw_boost` requires three paramters, `rs`, `b0` and `alpha`. The default values are in kpc. \n", + "- `compute_nfw_boost` requires two parameters to be specified, `rs` (non-optional) and `b0` (optional)\n", + "- `compute_powerlaw_boost` requires three paramters, `rs`(non-optional), `b0`(optional) and `alpha`(optional). \n", "\n", "Details on these boost models can be found [here](https://cluster-toolkit.readthedocs.io/en/latest/source/boostfactors.html)\n", "\n", @@ -246,23 +247,14 @@ "metadata": {}, "outputs": [], "source": [ - "nfw_boost = u.compute_nfw_boost(cl.DeltaSigma_profile[\"radius\"], rscale=1000, boost0=0.1)\n", + "r_proj = cl.DeltaSigma_profile[\"radius\"]\n", + "r_scale = 1000 # in kpc to match r_proj units\n", + "nfw_boost = u.compute_nfw_boost(r_proj, r_scale, boost0=0.1)\n", "\n", - "powerlaw_boost = u.compute_powerlaw_boost(\n", - " cl.DeltaSigma_profile[\"radius\"], rscale=1000, boost0=0.1, alpha=-1.0\n", + "powerlaw_boost = u.compute_powerlaw_boost(r_proj, r_scale, boost0=0.1, alpha=-1.0\n", ")" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "d10a993d-acbc-4525-807d-9f2851bf48d8", - "metadata": {}, - "outputs": [], - "source": [ - "print(u.compute_nfw_boost(1000, 1000, boost0=0.1))" - ] - }, { "cell_type": "markdown", "id": "2e356f74", @@ -278,11 +270,10 @@ "metadata": {}, "outputs": [], "source": [ - "plt.loglog(cl.DeltaSigma_profile[\"radius\"], nfw_boost, label=\"NFW boost factor\", marker='.')\n", - "plt.loglog(cl.DeltaSigma_profile[\"radius\"], powerlaw_boost, label=\"Powerlaw boost factor\")\n", + "plt.loglog(r_proj, nfw_boost, label=\"NFW boost factor\", marker='.')\n", + "plt.loglog(r_proj, powerlaw_boost, label=\"Powerlaw boost factor\", marker='.')\n", "plt.xlabel(\"Radius [kpc]\")\n", "plt.ylabel(\"$B(R)$\")\n", - "plt.axvline(1000)\n", "plt.legend()\n", "plt.show()" ] @@ -408,13 +399,13 @@ " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", " \"powerlaw_boost\", # boost_model\n", - " 1000, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", + " r_scale, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", ")\n", "Sigma_corrected_nfw_boost = u.correct_with_boost_model(\n", " cl.DeltaSigma_profile[\"radius\"],\n", " cl.DeltaSigma_profile[\"DeltaSigma_tan\"],\n", " \"nfw_boost\", # boost_model\n", - " 1000, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", + " r_scale, # boost_rscale (in units of cl.DeltaSigma_profile[\"radius\"])\n", ")\n", "\n", "plt.errorbar(\n", From a3d4ae0bbf1a21074293837db06112a0a6e5e161 Mon Sep 17 00:00:00 2001 From: m-aguena Date: Fri, 18 Oct 2024 17:40:07 +0200 Subject: [PATCH 10/10] Update version to 1.14.1 --- clmm/__init__.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/clmm/__init__.py b/clmm/__init__.py index 64fee9cf1..0b688448d 100644 --- a/clmm/__init__.py +++ b/clmm/__init__.py @@ -26,4 +26,4 @@ ) from . import support -__version__ = "1.14.0" +__version__ = "1.14.1"