gh-188: add docstrings to all functions and tidy docs

docs/conf.py

@@ -25,9 +25,11 @@
 extensions = [
     "matplotlib.sphinxext.plot_directive",
     "nbsphinx",
+    "sphinx.ext.autodoc",
     "sphinx.ext.intersphinx",
     "sphinx.ext.napoleon",
     "sphinx.ext.viewcode",
+    "sphinx_autodoc_typehints",
     "sphinxcontrib.katex",
 ]
 
@@ -75,6 +77,14 @@
 # -- napoleon ----------------------------------------------------------------
 
 napoleon_google_docstring = False
+napoleon_use_rtype = False
+
+
+# -- sphinx-autodoc-typehints ------------------------------------------------
+
+always_use_bars_union = True
+typehints_use_rtype = False
+typehints_defaults = "comma"
 
 
 # -- plot_directive ----------------------------------------------------------

glass/core/algorithm.py

@@ -14,17 +14,31 @@ def nnls(
     maxiter: int | None = None,
 ) -> npt.NDArray[np.float64]:
     """
-    Compute a non-negative least squares solution.
+    _summary_.
 
-    Implementation of the algorithm due to [1] as described in [2].
-
-    References
+    Parameters
     ----------
-    * [1] Lawson, C. L. and Hanson, R. J. (1995), Solving Least Squares
-        Problems. doi: 10.1137/1.9781611971217
-    * [2] Bro, R. and De Jong, S. (1997), A fast
-        non-negativity-constrained least squares algorithm. J.
-        Chemometrics, 11, 393-401.
+    a
+        _description_
+    b
+        _description_
+    tol
+        _description_
+    maxiter
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    Raises
+    ------
+    ValueError
+        _description_
+    ValueError
+        _description_
+    ValueError
+        _description_
 
     """
     a = np.asanyarray(a)

glass/core/array.py

@@ -12,7 +12,14 @@
 def broadcast_first(
     *arrays: npt.NDArray[np.float64],
 ) -> tuple[npt.NDArray[np.float64], ...]:
-    """Broadcast arrays, treating the first axis as common."""
+    """
+    _summary_.
+
+    Returns
+    -------
+        _description_
+
+    """
     arrays = tuple(np.moveaxis(a, 0, -1) if np.ndim(a) else a for a in arrays)
     arrays = np.broadcast_arrays(*arrays)
     return tuple(np.moveaxis(a, -1, 0) if np.ndim(a) else a for a in arrays)
@@ -28,28 +35,11 @@ def broadcast_leading_axes(
     typing.Unpack[tuple[npt.NDArray[np.float64], ...]],
 ]:
     """
-    Broadcast all but the last N axes.
-
-    Returns the shape of the broadcast dimensions, and all input arrays
-    with leading axes matching that shape.
-
-    Examples
-    --------
-    Broadcast all dimensions of ``a``, all except the last dimension of
-    ``b``, and all except the last two dimensions of ``c``.
-
-    >>> a = 0
-    >>> b = np.zeros((4, 10))
-    >>> c = np.zeros((3, 1, 5, 6))
-    >>> dims, a, b, c = broadcast_leading_axes((a, 0), (b, 1), (c, 2))
-    >>> dims
-    (3, 4)
-    >>> a.shape
-    (3, 4)
-    >>> b.shape
-    (3, 4, 10)
-    >>> c.shape
-    (3, 4, 5, 6)
+    _summary_.
+
+    Returns
+    -------
+        _description_
 
     """
     shapes, trails = [], []
@@ -72,7 +62,29 @@ def ndinterp(  # noqa: PLR0913
     right: float | None = None,
     period: float | None = None,
 ) -> npt.NDArray[np.float64]:
-    """Interpolate multi-dimensional array over axis."""
+    """
+    _summary_.
+
+    Parameters
+    ----------
+    xp
+        _description_
+    fp
+        _description_
+    axis
+        _description_
+    left
+        _description_
+    right
+        _description_
+    period
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    """
     return np.apply_along_axis(
         partial(np.interp, x, xp),
         axis,
@@ -91,7 +103,21 @@ def trapz_product(
     ],
     axis: int = -1,
 ) -> npt.NDArray[np.float64]:
-    """Trapezoidal rule for a product of functions."""
+    """
+    _summary_.
+
+    Parameters
+    ----------
+    f
+        _description_
+    axis
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    """
     x, _ = f
     for x_, _ in ff:
         x = np.union1d(
@@ -114,7 +140,25 @@ def cumtrapz(
     dtype: type | None = None,
     out: npt.NDArray[np.float64] | None = None,
 ) -> npt.NDArray[np.float64]:
-    """Cumulative trapezoidal rule along last axis."""
+    """
+    Cumulative trapezoidal rule along last axis.
+
+    Parameters
+    ----------
+    f
+        _description_
+    x
+        _description_
+    dtype
+        _description_
+    out
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    """
     if out is None:
         out = np.empty_like(f, dtype=dtype)
 

glass/ext/__init__.py

@@ -8,6 +8,21 @@
 
 
 def _extend_path(path: list[str], name: str) -> list[str]:
+    """
+    _summary_.
+
+    Parameters
+    ----------
+    path
+        _description_
+    name
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    """
     import os.path
     from pkgutil import extend_path
 

glass/fields.py

@@ -210,26 +210,31 @@ def generate_gaussian(
     rng: np.random.Generator | None = None,
 ) -> collections.abc.Generator[npt.NDArray[np.complex128]]:
     """
-    Sample Gaussian random fields from Cls iteratively.
+    _summary_.
 
-    A generator that iteratively samples HEALPix maps of Gaussian random fields
-    with the given angular power spectra ``gls`` and resolution parameter
-    ``nside``.
-
-    The optional argument ``ncorr`` can be used to artificially limit now many
-    realised fields are correlated. This saves memory, as only `ncorr` previous
-    fields need to be kept.
-
-    The ``gls`` array must contain the auto-correlation of each new field
-    followed by the cross-correlations with all previous fields in reverse
-    order::
-
-        gls = [gl_00,
-               gl_11, gl_10,
-               gl_22, gl_21, gl_20,
-               ...]
-
-    Missing entries can be set to ``None``.
+    Parameters
+    ----------
+    gls
+        _description_
+    nside
+        _description_
+    ncorr
+        _description_
+    rng
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    Yields
+    ------
+        _description_
+
+    Raises
+    ------
+    ValueError
+        _description_
 
     """
     # get the default RNG if not given
@@ -295,7 +300,31 @@ def generate_lognormal(
     ncorr: int | None = None,
     rng: np.random.Generator | None = None,
 ) -> collections.abc.Generator[npt.NDArray[np.complex128]]:
-    """Sample lognormal random fields from Gaussian Cls iteratively."""
+    """
+    _summary_.
+
+    Parameters
+    ----------
+    gls
+        _description_
+    nside
+        _description_
+    shift
+        _description_
+    ncorr
+        _description_
+    rng
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    Yields
+    ------
+        _description_
+
+    """
     for i, m in enumerate(generate_gaussian(gls, nside, ncorr=ncorr, rng=rng)):
         # compute the variance of the auto-correlation
         gl = gls[i * (i + 1) // 2]
@@ -323,21 +352,22 @@ def getcl(
     lmax: int | None = None,
 ) -> npt.NDArray[np.float64]:
     """
-    Return a specific angular power spectrum from an array.
-
-    Return the angular power spectrum for indices *i* and *j* from an
-    array in *GLASS* ordering.
+    _summary_.
 
     Parameters
     ----------
-    cls:
-        List of angular power spectra in *GLASS* ordering.
-    i:
-        Indices to return.
-    j:
-        Indices to return.
-    lmax:
-        Truncate the returned spectrum at this mode number.
+    cls
+        _description_
+    i
+        _description_
+    j
+        _description_
+    lmax
+        _description_
+
+    Returns
+    -------
+        _description_
 
     """
     if j > i:
@@ -358,28 +388,30 @@ def effective_cls(
     *,
     lmax: int | None = None,
 ) -> npt.NDArray[np.float64]:
-    r"""
-    Compute effective angular power spectra from weights.
-
-    Computes a linear combination of the angular power spectra *cls*
-    using the factors provided by *weights1* and *weights2*. Additional
-    axes in *weights1* and *weights2* produce arrays of spectra.
-
-    Returns a dictionary of effective angular power spectra, where keys
-    correspond to the leading axes of *weights1* and *weights2*.
+    """
+    _summary_.
 
     Parameters
     ----------
-    cls:
-        Angular matter power spectra to combine, in *GLASS* ordering.
-    weights1:
-        Weight factors for spectra. The first axis must be equal to
-        the number of fields.
-    weights2:
-        Second set of weights. If not given, *weights1* is used.
-    lmax:
-        Truncate the angular power spectra at this mode number. If not
-        given, the longest input in *cls* will be used.
+    cls
+        _description_
+    weights1
+        _description_
+    weights2
+        _description_
+    lmax
+        _description_
+
+    Returns
+    -------
+        _description_
+
+    Raises
+    ------
+    ValueError
+        _description_
+    ValueError
+        _description_
 
     """
     from itertools import combinations_with_replacement, product