Merge branch 'v1.0.0' of github.com:Cosmoglobe/zodipy into v1.0.0

mkdocs.yml

@@ -8,12 +8,12 @@ theme:
     - media: "(prefers-color-scheme: light)"
       scheme: default
       toggle:
-        icon: material/lightbulb
+        icon: material/weather-sunny
         name: Switch to dark mode
     - media: "(prefers-color-scheme: dark)"
       scheme: slate
       toggle:
-        icon:  material/lightbulb-outline
+        icon:  material/weather-night
         name: Switch to light mode
   icon:
     logo: material/book-open-page-variant
@@ -30,6 +30,7 @@ theme:
     - navigation.tabs
     - navigation.instant
     - navigation.top
+    - navigation.tabs.sticky
 markdown_extensions:
   - markdown_include.include:
       base_path: docs

tests/test_evaluate.py

@@ -1,7 +1,5 @@
 from __future__ import annotations
 
-import multiprocessing
-
 import numpy as np
 import pytest
 from astropy import coordinates as coords
@@ -52,7 +50,10 @@ def test_evaluate(
     obs: units.Quantity | str,
 ) -> None:
     """Test that the evaluate function works for valid user input."""
-    assert model.evaluate(sky_coord, obspos=obs).size == sky_coord.size
+    emission = model.evaluate(sky_coord, obspos=obs)
+    assert emission.size == sky_coord.size
+    assert isinstance(emission, units.Quantity)
+    assert emission.unit == units.MJy / units.sr
 
 
 def test_invalid_sky_coord() -> None:
@@ -145,23 +146,7 @@ def test_contains_duplicates() -> None:
 
 def test_multiprocessing_nproc() -> None:
     """Test that the multiprocessing works with n_proc > 1."""
-    model_multi = Model(x=20 * units.micron, n_proc=4)
-    model = Model(x=20 * units.micron, n_proc=1)
-
-    lon = np.random.randint(low=0, high=360, size=10000)
-    lat = np.random.randint(low=-90, high=90, size=10000)
-    skycoord = SkyCoord(
-        lon,
-        lat,
-        unit=units.deg,
-        obstime=TEST_TIME,
-    )
-    emission_multi = model_multi.evaluate(skycoord)
-    emission = model.evaluate(skycoord)
-    assert_array_equal(emission_multi, emission)
-
-    model_multi = Model(x=75 * units.micron, n_proc=4, name="rrm-experimental")
-    model = Model(x=75 * units.micron, n_proc=1, name="rrm-experimental")
+    model = Model(x=20 * units.micron)
 
     lon = np.random.randint(low=0, high=360, size=10000)
     lat = np.random.randint(low=-90, high=90, size=10000)
@@ -171,15 +156,12 @@ def test_multiprocessing_nproc() -> None:
         unit=units.deg,
         obstime=TEST_TIME,
     )
-    emission_multi = model_multi.evaluate(skycoord)
-    emission = model.evaluate(skycoord)
-
+    emission_multi = model.evaluate(skycoord, nprocesses=4)
+    emission = model.evaluate(skycoord, nprocesses=1)
     assert_array_equal(emission_multi, emission)
 
+    model = Model(x=75 * units.micron, name="rrm-experimental")
 
-def test_multiprocessing_pool() -> None:
-    """Test that the multiprocessing works with n_proc > 1."""
-    model = Model(x=20 * units.micron, n_proc=1)
     lon = np.random.randint(low=0, high=360, size=10000)
     lat = np.random.randint(low=-90, high=90, size=10000)
     skycoord = SkyCoord(
@@ -188,17 +170,7 @@ def test_multiprocessing_pool() -> None:
         unit=units.deg,
         obstime=TEST_TIME,
     )
-    emission = model.evaluate(skycoord)
-
-    try:
-        from pytest_cov.embed import cleanup_on_sigterm
-    except ImportError:
-        pass
-    else:
-        cleanup_on_sigterm()
-
-    with multiprocessing.Pool(processes=4) as pool:
-        model_multi = Model(x=20 * units.micron, pool=pool)
-        emission_multi = model_multi.evaluate(skycoord)
+    emission_multi = model.evaluate(skycoord, nprocesses=4)
+    emission = model.evaluate(skycoord, nprocesses=1)
 
     assert_array_equal(emission_multi, emission)

zodipy/__init__.py

@@ -1,4 +1,3 @@
-from zodipy import comps, source_params
 from zodipy.model_registry import model_registry
 from zodipy.number_density import grid_number_density
 from zodipy.zodipy import Model
@@ -7,6 +6,4 @@
     "Model",
     "grid_number_density",
     "model_registry",
-    "comps",
-    "source_params",
 )

zodipy/brightness.py

@@ -20,8 +20,10 @@
 
 def kelsall_brightness_at_step(
     r: npt.NDArray[np.float64],
-    start: np.float64,
+    start: npt.NDArray[np.float64],
     stop: npt.NDArray[np.float64],
+    quad_root_0: float,
+    quad_root_n: float,
     X_obs: npt.NDArray[np.float64],
     u_los: npt.NDArray[np.float64],
     bp_interpolation_table: npt.NDArray[np.float64],
@@ -36,8 +38,11 @@ def kelsall_brightness_at_step(
     solar_irradiance: np.float64,
 ) -> npt.NDArray[np.float64]:
     """Kelsall uses common line of sight grid from obs to 5.2 AU."""
-    # Convert the quadrature range from [-1, 1] to the true ecliptic positions
-    R_los = ((stop - start) / 2) * r + (stop + start) / 2
+    # Convert the quadrature range from [0, inf] to the true ecliptic positions
+    a, b = start, stop
+    i, j = quad_root_0, quad_root_n
+    R_los = ((b - a) / (j - i)) * r + (j * a - i * b) / (j - i)
+
     X_los = R_los * u_los
     X_helio = X_los + X_obs
     R_helio = np.sqrt(X_helio[0] ** 2 + X_helio[1] ** 2 + X_helio[2] ** 2)
@@ -59,6 +64,8 @@ def rrm_brightness_at_step(
     r: npt.NDArray[np.float64],
     start: npt.NDArray[np.float64],
     stop: npt.NDArray[np.float64],
+    quad_root_0: float,
+    quad_root_n: float,
     X_obs: npt.NDArray[np.float64],
     u_los: npt.NDArray[np.float64],
     bp_interpolation_table: npt.NDArray[np.float64],
@@ -68,8 +75,11 @@ def rrm_brightness_at_step(
     calibration: np.float64,
 ) -> npt.NDArray[np.float64]:
     """RRM is implented with component specific line-of-sight grids."""
-    # Convert the quadrature range from [-1, 1] to the true ecliptic positions
-    R_los = ((stop - start) / 2) * r + (stop + start) / 2
+    # Convert the quadrature range from [0, inf] to the true ecliptic positions
+    a, b = start, stop
+    i, j = quad_root_0, quad_root_n
+    R_los = ((b - a) / (j - i)) * r + (j * a - i * b) / (j - i)
+
     X_los = R_los * u_los
     X_helio = X_los + X_obs
     R_helio = np.sqrt(X_helio[0] ** 2 + X_helio[1] ** 2 + X_helio[2] ** 2)

zodipy/line_of_sight.py

@@ -52,13 +52,13 @@
 COMPONENT_CUTOFFS = {**DIRBE_CUTOFFS, **RRM_CUTOFFS}
 
 
-def integrate_gauss_legendre(
+def integrate_gauss_laguerre(
     fn: Callable[[float], npt.NDArray[np.float64]],
     points: npt.NDArray[np.float64],
     weights: npt.NDArray[np.float64],
 ) -> npt.NDArray[np.float64]:
-    """Integrate a function using Gauss-Legendre quadrature."""
-    return np.squeeze(sum(fn(x) * w for x, w in zip(points, weights)))
+    """Integrate a function using Gauss-Laguerre quadrature."""
+    return np.squeeze(sum(np.exp(x) * fn(x) * w for x, w in zip(points, weights)))
 
 
 def get_sphere_intersection(