Repository: jlvdb/yet_another_wizz (language: python)

Issue:
Unittests
update tests


Discussion / hints:

Fix patch:
diff --git a/CHANGELOG.rst b/CHANGELOG.rst
index b4d082b7..abfee101 100644
--- a/CHANGELOG.rst
+++ b/CHANGELOG.rst
@@ -2,6 +2,19 @@ Change log
 ==========
 
 
+Version 2.5.5
+-------------
+
+- Implemented a uniform API for configuration classes
+- Deprecated ``AutoBinningConfig`` and ``ManualBinningConfig`` and combined them
+  in new ``BinningConfig`` class.
+- Adhere to python standards in data model (double underscore methods).
+- Added missing type checks for data concatenation.
+- Improved the unit test coverage.
+- Corrected default values for ``.is_compatible()`` methods.
+- Corrected some errors in the documentation.
+
+
 Version 2.5.4
 -------------
 
diff --git a/src/yaw/__init__.py b/src/yaw/__init__.py
index cffa8fbf..2077d152 100644
--- a/src/yaw/__init__.py
+++ b/src/yaw/__init__.py
@@ -28,15 +28,14 @@
 from yaw.core.cosmology import Scale
 from yaw.core.math import global_covariance
 from yaw.correlation import CorrData, CorrFunc, autocorrelate, crosscorrelate
-from yaw.deprecated import (
-    CorrelationData,
-    CorrelationFunction,
-    HistogramData,
-    PairCountResult,
-)
 from yaw.randoms import UniformRandoms
 from yaw.redshifts import HistData, RedshiftData
 
+# isort: split
+from yaw.deprecated.correlation.corrfuncs import CorrelationData, CorrelationFunction
+from yaw.deprecated.correlation.paircounts import PairCountResult
+from yaw.deprecated.redshifts import HistogramData
+
 __all__ = [
     "NewCatalog",
     "Configuration",
@@ -56,4 +55,4 @@
     "HistogramData",
     "PairCountResult",
 ]
-__version__ = "2.5.4"
+__version__ = "2.5.5"
diff --git a/src/yaw/catalogs/catalog.py b/src/yaw/catalogs/catalog.py
index 32e313fd..f06a1a63 100644
--- a/src/yaw/catalogs/catalog.py
+++ b/src/yaw/catalogs/catalog.py
@@ -1,7 +1,7 @@
 from __future__ import annotations
 
 import logging
-from abc import abstractclassmethod, abstractmethod, abstractproperty
+from abc import abstractmethod
 from collections.abc import Iterator
 from typing import TYPE_CHECKING, Any, TypeVar
 
@@ -117,10 +117,10 @@ def from_file(
                 "centers)"
             )
 
-        cls._logger.info(f"reading catalog file '{filepath}'")
+        cls._logger.info("reading catalog file '%s'", filepath)
         data = apd.read_auto(filepath, columns=columns, ext=file_ext, **kwargs)
         if sparse is not None:
-            cls._logger.debug(f"sparse sampling data {sparse}x")
+            cls._logger.debug("sparse sampling data %ix", sparse)
             data = data[::sparse]
         return cls(
             data,
@@ -133,13 +133,14 @@ def from_file(
             progress=progress,
         )
 
-    @abstractclassmethod
+    @classmethod
+    @abstractmethod
     def from_cache(cls, cache_directory: str, progress: bool = False) -> BaseCatalog:
         """Restore the catalogue from its cache directory.
 
         Catalogs should be instantiated through the factory class, see
         :meth:`yaw.catalogs.NewCatalog.from_cache`."""
-        cls._logger.info(f"restoring from cache directory '{cache_directory}'")
+        cls._logger.info("restoring from cache directory '%s'", cache_directory)
 
     def __repr__(self) -> str:
         name = self.__class__.__name__
@@ -160,7 +161,8 @@ def __len__(self) -> int:
     def __getitem__(self, item: int) -> Any:
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def ids(self) -> list[int]:
         """Return a list of unique patch indices in the catalog."""
         pass
@@ -215,27 +217,32 @@ def pos(self) -> CoordSky:
         """
         return CoordSky(self.ra, self.dec)
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def ra(self) -> NDArray[np.float_]:
         """Get an array of the right ascension values in radians."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def dec(self) -> NDArray[np.float_]:
         """Get an array of the declination values in radians."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def redshifts(self) -> NDArray[np.float_] | None:
         """Get the redshifts as array or ``None`` if not available."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def weights(self) -> NDArray[np.float_]:
         """Get the object weights as array or ``None`` if not available."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def patch(self) -> NDArray[np.int_]:
         """Get the patch indices of each object as array."""
         pass
@@ -250,7 +257,8 @@ def get_max_redshift(self) -> float:
         """Get the maximum redshift or ``None`` if not available."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def total(self) -> float:
         """Get the sum of weights or the number of objects if weights are not
         available."""
@@ -260,7 +268,8 @@ def get_totals(self) -> NDArray[np.float_]:
         """Get an array of the sum of weights or number of objects in each
         patch."""
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def centers(self) -> CoordSky:
         """Get a vector of sky coordinates of the patch centers in radians.
 
@@ -269,7 +278,8 @@ def centers(self) -> CoordSky:
         """
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def radii(self) -> DistSky:
         """Get a vector of angular separations in radians that describe the
         patch sizes.
@@ -337,8 +347,11 @@ def correlate(
         n1 = long_num_format(len(self))
         n2 = long_num_format(len(self) if other is None else len(other))
         self._logger.debug(
-            f"correlating with {'' if binned else 'un'}binned catalog "
-            f"({n1}x{n2}) in {config.binning.zbin_num} redshift bins"
+            "correlating with %sbinned catalog (%sx%s) in %d redshift bins",
+            "" if binned else "un",
+            n1,
+            n2,
+            config.binning.zbin_num,
         )
 
     @abstractmethod
diff --git a/src/yaw/catalogs/linkage.py b/src/yaw/catalogs/linkage.py
index a4a92df3..3042fa90 100644
--- a/src/yaw/catalogs/linkage.py
+++ b/src/yaw/catalogs/linkage.py
@@ -76,7 +76,7 @@ def from_setup(cls, config: Configuration, catalog: BaseCatalog) -> PatchLinkage
             max_query_radius = 0.0  # only relevant for cross-patch
         max_query_radius = DistSky(max_query_radius)
 
-        logger.debug(f"computing patch linkage with {max_query_radius=:.3e}")
+        logger.debug("computing patch linkage with %.3e", max_query_radius.values)
         centers_3d = catalog.centers.to_3d().values
         radii = catalog.radii.values
         # compute distance between all patch centers
@@ -90,7 +90,7 @@ def from_setup(cls, config: Configuration, catalog: BaseCatalog) -> PatchLinkage
         for id1, overlap in enumerate(overlaps):
             patch_pairs.extend((id1, id2) for id2 in np.where(overlap)[0])
         logger.debug(
-            f"found {len(patch_pairs)} patch links " f"for {catalog.n_patches} patches"
+            "found %d patch links for %d patches", len(patch_pairs), catalog.n_patches
         )
         return cls(patch_pairs)
 
diff --git a/src/yaw/catalogs/scipy/catalog.py b/src/yaw/catalogs/scipy/catalog.py
index 9377c209..a71e5844 100644
--- a/src/yaw/catalogs/scipy/catalog.py
+++ b/src/yaw/catalogs/scipy/catalog.py
@@ -165,7 +165,7 @@ def __init__(
                 raise FileNotFoundError(
                     f"patch directory does not exist: '{cache_directory}'"
                 )
-            self._logger.debug(f"using cache directory '{cache_directory}'")
+            self._logger.debug("using cache directory '%s'", cache_directory)
 
         # create new patches
         if patch_mode != "dividing":
@@ -176,21 +176,21 @@ def __init__(
                 patch_centers, patch_ids = create_patches(
                     position=position, n_patches=n_patches
                 )
-                log_msg = f"creating {n_patches} patches"
+                log_msg = "creating %i patches"
             else:
                 if isinstance(patch_centers, BaseCatalog):
                     patch_centers = patch_centers.centers.to_3d()
                 patch_ids = assign_patches(centers=patch_centers, position=position)
                 n_patches = len(patch_centers)
-                log_msg = f"applying {n_patches} patches from external data"
+                log_msg = "applying %i patches from external data"
             patch_name = "patch"  # the default name
             data[patch_name] = patch_ids
             centers = {pid: pos for pid, pos in enumerate(patch_centers)}
         else:
             n_patches = len(data[patch_name].unique())
-            log_msg = f"dividing data into {n_patches} predefined patches"
+            log_msg = "dividing data into %i predefined patches"
             centers = dict()  # this can be empty
-        self._logger.debug(log_msg)
+        self._logger.debug(log_msg, n_patches)
 
         # run groupby first to avoid any intermediate copies of full data
         n_obj_str = long_num_format(len(data))
diff --git a/src/yaw/catalogs/treecorr/catalog.py b/src/yaw/catalogs/treecorr/catalog.py
index 85509bdc..4aa4bcb6 100644
--- a/src/yaw/catalogs/treecorr/catalog.py
+++ b/src/yaw/catalogs/treecorr/catalog.py
@@ -128,7 +128,7 @@ def __init__(
                 raise FileNotFoundError(
                     f"patch directory does not exist: '{cache_directory}'"
                 )
-            self._logger.info(f"using cache directory '{cache_directory}'")
+            self._logger.info("using cache directory '%s'", cache_directory)
 
         if n_patches is not None:
             kwargs["npatch"] = n_patches
@@ -366,7 +366,7 @@ def correlate(
 
         # iterate the bins and compute the correlation
         self._logger.debug(
-            f"running treecorr on {config.backend.get_threads()} threads"
+            "running treecorr on %i threads", config.backend.get_threads()
         )
         for i, ((intv, bincat1), (_, bincat2)) in enumerate(zip(cats1, cats2)):
             if progress:
diff --git a/src/yaw/config/__init__.py b/src/yaw/config/__init__.py
index c0b80dd0..96cc2c8c 100644
--- a/src/yaw/config/__init__.py
+++ b/src/yaw/config/__init__.py
@@ -29,17 +29,21 @@
 
 # isort: split
 from yaw.config.backend import BackendConfig
-from yaw.config.binning import AutoBinningConfig, ManualBinningConfig
+from yaw.config.binning import BinningConfig
 from yaw.config.config import Configuration
 from yaw.config.resampling import ResamplingConfig
 from yaw.config.scales import ScalesConfig
 
+# isort: split
+from yaw.deprecated.config.binning import AutoBinningConfig, ManualBinningConfig
+
 __all__ = [
     "DEFAULT",
-    "BackendConfig",
     "AutoBinningConfig",
-    "ManualBinningConfig",
+    "BackendConfig",
+    "BinningConfig",
     "Configuration",
+    "ManualBinningConfig",
     "OPTIONS",
     "ResamplingConfig",
     "ScalesConfig",
diff --git a/src/yaw/config/abc.py b/src/yaw/config/abc.py
new file mode 100644
index 00000000..846dc89e
--- /dev/null
+++ b/src/yaw/config/abc.py
@@ -0,0 +1,53 @@
+from __future__ import annotations
+
+from abc import abstractmethod
+from typing import Any, Type, TypeVar
+
+from yaw.config import DEFAULT
+from yaw.core.abc import DictRepresentation
+
+T = TypeVar("T", bound="BaseConfig")
+
+
+class BaseConfig(DictRepresentation):
+    @classmethod
+    def create(cls: Type[T], **kwargs: Any) -> T:
+        """Create a new configuration object.
+
+        By default this is an alias for :meth:`__init__`. Configuration classes
+        that are hierarchical (i.e. contain configuration objects as attributes)
+        implement this method to provide a single constructer for its own and
+        its subclasses parameters.
+        """
+        return cls(**kwargs)
+
+    @abstractmethod
+    def modify(self: T, **kwargs: Any | DEFAULT.NotSet) -> T:
+        """Create a copy of the current configuration with updated parameter
+        values.
+
+        The method arguments are identical to :meth:`create`. Values that should
+        not be modified are by default represented by the special value
+        :obj:`~yaw.config.default.NotSet`.
+        """
+        conf_dict = self.to_dict()
+        conf_dict.update(
+            {key: value for key, value in kwargs.items() if value is not DEFAULT.NotSet}
+        )
+        return self.__class__.from_dict(conf_dict)
+
+    @classmethod
+    def from_dict(
+        cls: Type[T],
+        the_dict: dict[str, Any],
+        **kwargs: dict[str, Any],  # passing additional constructor data
+    ) -> T:
+        """Create a class instance from a dictionary representation of the
+        minimally required data.
+
+        Args:
+            the_dict (:obj:`dict`):
+                Dictionary containing the data.
+            **kwargs: Additional data needed to construct the class instance.
+        """
+        return cls.create(**the_dict)
diff --git a/src/yaw/config/backend.py b/src/yaw/config/backend.py
index 99b620a5..a65067f2 100644
--- a/src/yaw/config/backend.py
+++ b/src/yaw/config/backend.py
@@ -3,15 +3,15 @@
 import os
 from dataclasses import dataclass, field
 
-from yaw.config import default as DEFAULT
-from yaw.core.abc import DictRepresentation
+from yaw.config import DEFAULT
+from yaw.config.abc import BaseConfig
 from yaw.core.docs import Parameter
 
 __all__ = ["BackendConfig"]
 
 
 @dataclass(frozen=True)
-class BackendConfig(DictRepresentation):
+class BackendConfig(BaseConfig):
     """Configuration of backends used for correlation measurements.
 
     Args:
@@ -60,6 +60,16 @@ def __post_init__(self) -> None:
         if self.thread_num is None:
             object.__setattr__(self, "thread_num", os.cpu_count())
 
+    def modify(
+        self,
+        thread_num: int | None = DEFAULT.NotSet,
+        crosspatch: bool = DEFAULT.NotSet,
+        rbin_slop: float = DEFAULT.NotSet,
+    ) -> BackendConfig:
+        return super().modify(
+            thread_num=thread_num, crosspatch=crosspatch, rbin_slop=rbin_slop
+        )
+
     def get_threads(self, max=None) -> int:
         """Get the number of threads for parallel processing.
 
diff --git a/src/yaw/config/binning.py b/src/yaw/config/binning.py
index 675fa05e..c572dba2 100644
--- a/src/yaw/config/binning.py
+++ b/src/yaw/config/binning.py
@@ -6,10 +6,9 @@
 
 import numpy as np
 
-from yaw.config import OPTIONS
-from yaw.config import default as DEFAULT
+from yaw.config import DEFAULT, OPTIONS
+from yaw.config.abc import BaseConfig
 from yaw.config.utils import ConfigError
-from yaw.core.abc import DictRepresentation
 from yaw.core.cosmology import BinFactory, TypeCosmology
 from yaw.core.docs import Parameter
 from yaw.core.math import array_equal
@@ -17,34 +16,12 @@
 if TYPE_CHECKING:  # pragma: no cover
     from numpy.typing import NDArray
 
-__all__ = ["AutoBinningConfig", "ManualBinningConfig", "make_binning_config"]
+__all__ = ["BinningConfig"]
 
 
-class BaseBinningConfig(DictRepresentation):
-    """Base class for redshift binning configuration."""
-
-    zbins: NDArray[np.float_]
-    """Edges of redshift bins."""
-    method: str
-    """Method used to create redshift binning, ``manual`` or either of
-    :obj:`~yaw.config.options.Options.binning`."""
-
-    def __repr__(self) -> str:
-        name = self.__class__.__name__
-        zbin_num = self.zbin_num
-        z = f"{self.zmin:.3f}...{self.zmax:.3f}"
-        method = self.method
-        return f"{name}({zbin_num=}, {z=}, {method=})"
-
-
-@dataclass(frozen=True, repr=False)
-class ManualBinningConfig(BaseBinningConfig):
-    """Configuration that specifies a manual redshift binning.
-
-    Args:
-        zbins (:obj:`NDArray`):
-            Edges of redshift bins, must increase monotonically.
-    """
+@dataclass(frozen=True)
+class BinningConfig(BaseConfig):
+    """TODO"""
 
     zbins: NDArray[np.float_] = field(
         metadata=Parameter(
@@ -53,73 +30,17 @@ class ManualBinningConfig(BaseBinningConfig):
             help="list of custom redshift bin edges, if method is set to 'manual'",
         )
     )
-
-    def __post_init__(self) -> None:
-        if len(self.zbins) < 2:
-            raise ConfigError("'zbins' must have at least two edges")
-        BinFactory.check(self.zbins)
-        object.__setattr__(self, "zbins", np.asarray(self.zbins))
-
-    def __eq__(self, other: ManualBinningConfig) -> bool:
-        if not array_equal(self.zbins, other.zbins):
-            return False
-        return True
-
-    @property
-    def method(self) -> str:
-        """Method used to create redshift binning, always ``manual``."""
-        return "manual"
-
-    @property
-    def zmin(self) -> float:
-        """Lowest redshift bin edge."""
-        return float(self.zbins[0])
-
-    @property
-    def zmax(self) -> float:
-        """Highest redshift bin edge."""
-        return float(self.zbins[-1])
-
-    @property
-    def zbin_num(self) -> int:
-        """Number of redshift bins."""
-        return len(self.zbins) - 1
-
-    @classmethod
-    def from_dict(cls, the_dict: dict[str, Any], **kwargs) -> ManualBinningConfig:
-        new_dict = {key: val for key, val in the_dict.items()}
-        zbins = new_dict.pop("zbins")
-        return cls(np.asarray(zbins), **new_dict)
-
-    def to_dict(self) -> dict[str, Any]:
-        return dict(method=self.method, zbins=self.zbins.tolist())
-
-
-@dataclass(frozen=True, repr=False)
-class AutoBinningConfig(BaseBinningConfig):
-    """Configuration that generates a redshift binning.
-
-    To generate a redshift binning use the :meth:`generate` method.
-
-    Args:
-        zbins (:obj:`NDArray`):
-            Edges of redshift bins, must increase monotonically.
-        method (:obj:`str`):
-            Method used to create redshift binning, either of
-            :obj:`~yaw.config.options.Options.binning`.
-    """
-
-    zbins: NDArray[np.float_]
+    """Edges of redshift bins."""
     method: str = field(
-        default=DEFAULT.Configuration.binning.method,
+        default=DEFAULT.Binning.method,
         metadata=Parameter(
             choices=OPTIONS.binning,
             help="redshift binning method, 'logspace' means equal size in log(1+z)",
             default_text="(default: %(default)s)",
         ),
     )
-    """Method used to create redshift binning, see
-    :obj:`~yaw.config.options.Options.binning`."""
+    """Method used to create redshift binning (defaults to 'manual' if `zbins`
+    are provided), see :obj:`~yaw.config.options.Options.binning` for options."""
     zmin: float = field(
         init=False,
         metadata=Parameter(
@@ -139,7 +60,7 @@ class AutoBinningConfig(BaseBinningConfig):
     )
     """Highest redshift bin edge."""
     zbin_num: int = field(
-        default=DEFAULT.AutoBinning.zbin_num,
+        default=DEFAULT.Binning.zbin_num,
         init=False,
         metadata=Parameter(
             type=int,
@@ -150,26 +71,60 @@ class AutoBinningConfig(BaseBinningConfig):
     """Number of redshift bins."""
 
     def __post_init__(self) -> None:
+        if len(self.zbins) < 2:
+            raise ConfigError("'zbins' must have at least two edges")
+        object.__setattr__(self, "zbins", np.asarray(self.zbins))
+        BinFactory.check(self.zbins)
         object.__setattr__(self, "zmin", float(self.zbins[0]))
         object.__setattr__(self, "zmax", float(self.zbins[-1]))
         object.__setattr__(self, "zbin_num", len(self.zbins) - 1)
 
+    def __eq__(self, other) -> bool:
+        if not isinstance(other, self.__class__):
+            return False
+        if self.method != other.method or not array_equal(self.zbins, other.zbins):
+            return False
+        return True
+
+    def __repr__(self) -> str:
+        name = self.__class__.__name__
+        zbin_num = self.zbin_num
+        z = f"{self.zmin:.3f}...{self.zmax:.3f}"
+        method = self.method
+        return f"{name}({zbin_num=}, {z=}, {method=})"
+
+    @property
+    def is_manual(self) -> bool:
+        return self.method == "manual"
+
     @classmethod
-    def generate(
+    def create(
         cls,
-        zmin: float,
-        zmax: float,
-        zbin_num: int = DEFAULT.AutoBinning.zbin_num,
-        method: str = DEFAULT.AutoBinning.method,
-        cosmology: TypeCosmology | None = None,
-    ) -> AutoBinningConfig:
-        """Generate a new redshift binning configuration.
-
-        Generates a specified number of bins between a lower and upper redshift
-        limit. The spacing of the bins depends the generation method, the
-        default is a linear spacing.
+        *,
+        zbins: NDArray[np.float_] | None = None,
+        zmin: float | None = None,
+        zmax: float | None = None,
+        zbin_num: int = DEFAULT.Binning.zbin_num,
+        method: str = DEFAULT.Binning.method,
+        cosmology: TypeCosmology | str | None = None,
+    ) -> BinningConfig:
+        """Create a new redshift binning configuration.
+
+        If redshift bins (``zbins``) are provided, ``method`` is set to
+        ``"manual"`` and the bins edges are used. If ``zmin`` and ``zmax`` are
+        provided, instead generates a specified number of bins between this
+        lower and upper redshift limit. The spacing of the bins depends the
+        generation method, the default is a linear spacing.
+
+        .. Note::
+
+            The ``cosmology`` parameter is only used when generating binnings
+            that require cosmological distance computations.
 
         Args:
+            zbins (:obj:`NDArray[np.float_]`):
+                Monotonically increasing redshift bin edges, including the upper
+                edge (ignored if ``zmin`` and ``zmax`` are provided).
             zmin (:obj:`float`):
                 Minimum redshift, lowest redshift edge.
             zmax (:obj:`float`):
@@ -185,24 +140,75 @@ def generate(
                 model, refer to :mod:`yaw.core.cosmology`.
 
         Returns:
-            :obj:`AutoBinningConfig`
+            :obj:`BinningConfig`
+        """
+        auto_args_set = (zmin is not None, zmax is not None)
+        manual_args_set = (zbins is not None,)
+        if not all(manual_args_set) and not all(auto_args_set):
+            raise ConfigError("either 'zbins' or 'zmin' and 'zmax' are required")
+
+        elif all(auto_args_set):  # generate zbins
+            if all(manual_args_set):
+                warnings.warn(
+                    "'zbins' set but ignored since 'zmin' and 'zmax' are provided"
+                )
+            if not isinstance(method, str):
+                raise ValueError("'method' must a string")
+            zbins = BinFactory(zmin, zmax, zbin_num, cosmology).get(method)
+
+        else:  # use provided zbins
+            method = "manual"
+
+        return cls(zbins=zbins, method=method)
+
+    def modify(
+        self,
+        zbins: NDArray[np.float_] | None = DEFAULT.NotSet,
+        zmin: float | None = DEFAULT.NotSet,
+        zmax: float | None = DEFAULT.NotSet,
+        zbin_num: int = DEFAULT.NotSet,
+        method: str = DEFAULT.NotSet,
+        cosmology: TypeCosmology | str | None = None,
+    ) -> BinningConfig:
+        """Create a copy of the current configuration with updated parameter
+        values.
+
+        The method arguments are identical to :meth:`create`. Values that should
+        not be modified are by default represented by the special value
+        :obj:`~yaw.config.default.NotSet`.
+
+        .. Note::
+
+            The ``cosmology`` parameter is only used when generating binnings
+            that require cosmological distance computations.
         """
-        if not isinstance(method, str):
-            raise ValueError("'method' must a string")
-        zbins = BinFactory(zmin, zmax, zbin_num, cosmology).get(method)
-        return cls(zbins, method)
+        if zbins is not DEFAULT.NotSet:
+            kwargs = dict(zbins=zbins)
+        else:
+            inputs = dict(zmin=zmin, zmax=zmax, zbin_num=zbin_num, method=method)
+            mods = {k: v for k, v in inputs.items() if v is not DEFAULT.NotSet}
+            if self.is_manual:  # use every input as parameter
+                kwargs = mods
+            else:  # keep the existing parameters and update with inputs
+                kwargs = self.to_dict()
+                kwargs.update(mods)
+        return self.__class__.create(**kwargs, cosmology=cosmology)
 
-    def __eq__(self, other: ManualBinningConfig) -> bool:
-        if self.method != other.method:
-            return False
-        if not array_equal(self.zbins, other.zbins):
-            return False
-        return True
+    def to_dict(self) -> dict[str, Any]:
+        if self.is_manual:
+            return dict(method=self.method, zbins=self.zbins.tolist())
+        else:
+            return dict(
+                zmin=self.zmin,
+                zmax=self.zmax,
+                zbin_num=self.zbin_num,
+                method=self.method,
+            )
 
     @classmethod
     def from_dict(
         cls, the_dict: dict[str, Any], cosmology: TypeCosmology | None = None
-    ) -> AutoBinningConfig:
+    ) -> BinningConfig:
         """Create a class instance from a dictionary representation of the
         minimally required data.
 
@@ -213,42 +219,4 @@ def from_dict(
                 Cosmological model used for distance calculations. For a custom
                 model, refer to :mod:`yaw.core.cosmology`.
         """
-        return cls.generate(**the_dict, cosmology=cosmology)
-
-    def to_dict(self) -> dict[str, Any]:
-        return dict(
-            zmin=self.zmin, zmax=self.zmax, zbin_num=self.zbin_num, method=self.method
-        )
-
-
-def make_binning_config(
-    *,
-    cosmology: TypeCosmology | str | None = None,
-    zmin: float | None = None,
-    zmax: float | None = None,
-    zbin_num: int = DEFAULT.AutoBinning.zbin_num,
-    method: str = DEFAULT.AutoBinning.method,
-    zbins: NDArray[np.float_] | None = None,
-) -> ManualBinningConfig | AutoBinningConfig:
-    """
-    Helper function to construct a binning configuration.
-
-    The ``cosmology`` argument is always required. If redshift bins (``zbins``)
-    are provided, a :obj:`ManualBinningConfig` is returned, otherwise an
-    :obj:`AutoBinningConfig`. Issues a warning if any argument is set but
-    ignored by the returned configuration instance.
-    """
-    auto_args_set = (zmin is not None, zmax is not None)
-    manual_args_set = (zbins is not None,)
-    if not all(manual_args_set) and not all(auto_args_set):
-        raise ConfigError("either 'zbins' or 'zmin' and 'zmax' are required")
-    elif all(auto_args_set):
-        if all(manual_args_set):
-            warnings.warn(
-                "'zbins' set but ignored since 'zmin' and 'zmax' are provided"
-            )
-        return AutoBinningConfig.generate(
-            zmin=zmin, zmax=zmax, zbin_num=zbin_num, method=method, cosmology=cosmology
-        )
-    else:
-        return ManualBinningConfig(zbins)
+        return cls.create(**the_dict, cosmology=cosmology)
diff --git a/src/yaw/config/config.py b/src/yaw/config/config.py
index 408f62df..32866031 100644
--- a/src/yaw/config/config.py
+++ b/src/yaw/config/config.py
@@ -11,14 +11,10 @@
 from yaw.config import OPTIONS
 from yaw.config import default as DEFAULT
 from yaw.config import utils
+from yaw.config.abc import BaseConfig
 from yaw.config.backend import BackendConfig
-from yaw.config.binning import (
-    AutoBinningConfig,
-    ManualBinningConfig,
-    make_binning_config,
-)
+from yaw.config.binning import BinningConfig
 from yaw.config.scales import ScalesConfig
-from yaw.core.abc import DictRepresentation
 from yaw.core.cosmology import TypeCosmology, get_default_cosmology, r_kpc_to_angle
 from yaw.core.docs import Parameter
 
@@ -35,26 +31,8 @@
 logger = logging.getLogger(__name__)
 
 
-def update_if_set(the_dict: dict, key: str, value: Any | DEFAULT.NotSet) -> None:
-    if value is not DEFAULT.NotSet:
-        the_dict[key] = value
-
-
-def update_auto_binning_if_set(
-    the_dict: dict,
-    zmin: float | None = DEFAULT.NotSet,
-    zmax: float | None = DEFAULT.NotSet,
-    zbin_num: int | None = DEFAULT.NotSet,
-    method: str | None = DEFAULT.NotSet,
-) -> None:
-    update_if_set(the_dict, "zmin", zmin)
-    update_if_set(the_dict, "zmax", zmax)
-    update_if_set(the_dict, "zbin_num", zbin_num)
-    update_if_set(the_dict, "method", method)
-
-
 @dataclass(frozen=True)
-class Configuration(DictRepresentation):
+class Configuration(BaseConfig):
     """The central configration for correlation measurements.
 
     Bundles the configuration of measurement scales, redshift binning, and
@@ -92,7 +70,7 @@ class Configuration(DictRepresentation):
 
     scales: ScalesConfig
     """The configuration of the measurement scales."""
-    binning: AutoBinningConfig | ManualBinningConfig
+    binning: BinningConfig
     """The redshift binning configuration."""
     backend: BackendConfig = field(default_factory=BackendConfig)
     """The backend-specific configuration."""
@@ -196,8 +174,10 @@ def create(
             :obj:`Configuration`
         """
         cosmology = utils.parse_cosmology(cosmology)
-        scales = ScalesConfig(rmin=rmin, rmax=rmax, rweight=rweight, rbin_num=rbin_num)
-        binning = make_binning_config(
+        scales = ScalesConfig.create(
+            rmin=rmin, rmax=rmax, rweight=rweight, rbin_num=rbin_num
+        )
+        binning = BinningConfig.create(
             cosmology=cosmology,
             zmin=zmin,
             zmax=zmax,
@@ -205,7 +185,7 @@ def create(
             method=method,
             zbins=zbins,
         )
-        backend = BackendConfig(
+        backend = BackendConfig.create(
             thread_num=thread_num, crosspatch=crosspatch, rbin_slop=rbin_slop
         )
         return cls(scales=scales, binning=binning, backend=backend, cosmology=cosmology)
@@ -230,50 +210,27 @@ def modify(
         crosspatch: bool | None = DEFAULT.NotSet,
         rbin_slop: float | None = DEFAULT.NotSet,
     ) -> Configuration:
-        """Create a copy of the current configuration with updated parameter
-        values.
-
-        The method arguments are identical to :meth:`create`. Values that should
-        not be modified are by default represented by the special value
-        :obj:`~yaw.config.default.NotSet`.
-        """
-        config = self.to_dict()
-
-        # cosmology
-        if cosmology is not DEFAULT.NotSet:
-            if isinstance(cosmology, str):
-                cosmology = utils.yaml_to_cosmology(cosmology)
-            config["cosmology"] = utils.cosmology_to_yaml(cosmology)
-
-        # ScalesConfig
-        update_if_set(config["scales"], "rmin", rmin)
-        update_if_set(config["scales"], "rmax", rmax)
-        update_if_set(config["scales"], "rweight", rweight)
-        update_if_set(config["scales"], "rbin_num", rbin_num)
-
-        # AutoBinningConfig / ManualBinningConfig
-        if isinstance(self.binning, AutoBinningConfig):
-            if zbins is not DEFAULT.NotSet:
-                config["binning"] = dict(zbins=zbins)  # remove other params
-            else:
-                update_auto_binning_if_set(
-                    config["binning"], zmin, zmax, zbin_num, method
-                )
-        else:  # is ManualBinningConfig
-            if zbins is not DEFAULT.NotSet:
-                config["binning"]["zbins"] = zbins
-            else:
-                config["binning"] = dict()  # clear the binning parameters
-                update_auto_binning_if_set(
-                    config["binning"], zmin, zmax, zbin_num, method
-                )
-
-        # BackendConfig
-        update_if_set(config["backend"], "thread_num", thread_num)
-        update_if_set(config["backend"], "crosspatch", crosspatch)
-        update_if_set(config["backend"], "rbin_slop", rbin_slop)
-
-        return self.__class__.from_dict(config)
+        if cosmology is DEFAULT.NotSet:
+            cosmology = self.cosmology
+        elif isinstance(cosmology, str):
+            cosmology = utils.yaml_to_cosmology(cosmology)
+        scales = self.scales.modify(
+            rmin=rmin, rmax=rmax, rweight=rweight, rbin_num=rbin_num
+        )
+        binning = self.binning.modify(
+            zmin=zmin,
+            zmax=zmax,
+            method=method,
+            zbin_num=zbin_num,
+            zbins=zbins,
+            cosmology=cosmology,
+        )
+        backend = self.backend.modify(
+            thread_num=thread_num, crosspatch=crosspatch, rbin_slop=rbin_slop
+        )
+        return self.__class__(
+            cosmology=cosmology, scales=scales, binning=binning, backend=backend
+        )
 
     @deprecated(reason="no longer maintained", version="2.5.3")
     def plot_scales(
@@ -332,17 +289,19 @@ def from_dict(cls, the_dict: dict[str, Any], **kwargs) -> Configuration:
         )
         # parse the required subgroups
         try:
-            scales = ScalesConfig.from_dict(config.pop("scales"))
+            scales_dict = config.pop("scales")
+            scales = ScalesConfig.from_dict(scales_dict)
         except (TypeError, KeyError) as e:
             utils.parse_section_error(e, "scales")
         try:
-            binning_conf = config.pop("binning")
-            binning = make_binning_config(cosmology=cosmology, **binning_conf)
+            binning_dict = config.pop("binning")
+            binning = BinningConfig.from_dict(binning_dict, cosmology=cosmology)
         except (TypeError, KeyError) as e:
             utils.parse_section_error(e, "binning")
         # parse the optional subgroups
         try:
-            backend = BackendConfig.from_dict(config.pop("backend"))
+            backend_dict = config.pop("backend")
+            backend = BackendConfig.from_dict(backend_dict)
         except KeyError:
             backend = BackendConfig()
         except TypeError as e:
@@ -374,7 +333,7 @@ def from_yaml(cls, path: TypePathStr) -> Configuration:
         Returns:
             :obj:`Configuration`
         """
-        logger.info(f"reading configuration file '{path}'")
+        logger.info("reading configuration file '%s'", path)
         with open(str(path)) as f:
             config = yaml.safe_load(f.read())
         return cls.from_dict(config)
@@ -386,7 +345,7 @@ def to_yaml(self, path: TypePathStr) -> None:
             path (:obj:`pathlib.Path`, :obj:`str`):
                 Path to which the YAML file is written.
         """
-        logger.info(f"writing configuration file '{path}'")
+        logger.info("writing configuration file '%s'", path)
         string = yaml.dump(self.to_dict())
         with open(str(path), "w") as f:
             f.write(string)
diff --git a/src/yaw/config/default.py b/src/yaw/config/default.py
index 8431bff4..0d1abb6f 100644
--- a/src/yaw/config/default.py
+++ b/src/yaw/config/default.py
@@ -6,7 +6,7 @@
 __all__ = [
     "NotSet",
     "Scales",
-    "AutoBinning",
+    "Binning",
     "Backend",
     "Configuration",
     "Resampling",
@@ -31,11 +31,14 @@ class Scales:
     rbin_num = 50
 
 
-class AutoBinning:
+class Binning:
     method = "linear"
     zbin_num = 30
 
 
+AutoBinning = Binning  # keep for backwards compatibilty
+
+
 class Backend:
     thread_num = None
     crosspatch = True
@@ -44,7 +47,7 @@ class Backend:
 
 class Configuration:
     scales = Scales
-    binning = AutoBinning
+    binning = Binning
     backend = Backend
     cosmology = "Planck15"
 
diff --git a/src/yaw/config/resampling.py b/src/yaw/config/resampling.py
index e4a1e55a..b88fb4f8 100644
--- a/src/yaw/config/resampling.py
+++ b/src/yaw/config/resampling.py
@@ -7,8 +7,8 @@
 
 from yaw.config import OPTIONS
 from yaw.config import default as DEFAULT
+from yaw.config.abc import BaseConfig
 from yaw.config.utils import ConfigError
-from yaw.core.abc import DictRepresentation
 
 if TYPE_CHECKING:  # pragma: no cover
     from numpy.typing import NDArray
@@ -17,7 +17,7 @@
 
 
 @dataclass(frozen=True)
-class ResamplingConfig(DictRepresentation):
+class ResamplingConfig(BaseConfig):
     """Configuration for error estimation from spatial resampling.
 
     Used for all functions and methods that use spatial patches for error
@@ -62,6 +62,22 @@ def __post_init__(self) -> None:
                 f"invalid resampling method '{self.method}', must either of {opts}"
             )
 
+    def modify(
+        self,
+        method: str = DEFAULT.NotSet,
+        crosspatch: bool = DEFAULT.NotSet,
+        n_boot: int = DEFAULT.NotSet,
+        global_norm: bool = DEFAULT.NotSet,
+        seed: int = DEFAULT.NotSet,
+    ) -> ResamplingConfig:
+        return super().modify(
+            method=method,
+            crosspatch=crosspatch,
+            n_boot=n_boot,
+            global_norm=global_norm,
+            seed=seed,
+        )
+
     @property
     def n_patches(self) -> int | None:
         """The number of spatial patches for which this configuratin is valid.
diff --git a/src/yaw/config/scales.py b/src/yaw/config/scales.py
index da12dab8..112fc9c6 100644
--- a/src/yaw/config/scales.py
+++ b/src/yaw/config/scales.py
@@ -8,8 +8,8 @@
 from deprecated import deprecated
 
 from yaw.config import default as DEFAULT
+from yaw.config.abc import BaseConfig
 from yaw.config.utils import ConfigError
-from yaw.core.abc import DictRepresentation
 from yaw.core.cosmology import Scale
 from yaw.core.docs import Parameter
 from yaw.core.math import array_equal
@@ -21,7 +21,7 @@
 
 
 @dataclass(frozen=True)
-class ScalesConfig(DictRepresentation):
+class ScalesConfig(BaseConfig):
     """Configuration of scales used for correlation measurements.
 
     Correlation functions are measured on one or many intervals
@@ -120,14 +120,14 @@ def __post_init__(self) -> None:
         else:
             raise ConfigError("'rmin' and 'rmax' must be both sequences or float")
 
-    def __eq__(self, other: ScalesConfig) -> bool:
-        if not array_equal(self.as_array(), other.as_array()):
-            return False
-        if self.rweight != other.rweight:
-            return False
-        if self.rbin_num != other.rbin_num:
-            return False
-        return True
+    def __eq__(self, other: object) -> bool:
+        if isinstance(other, self.__class__):
+            return (
+                array_equal(self.as_array(), other.as_array())
+                and self.rweight == other.rweight
+                and self.rbin_num == other.rbin_num
+            )
+        return NotImplemented
 
     def __getitem__(self, idx: int) -> Scale:
         scales = self.as_array()
@@ -137,6 +137,15 @@ def __iter__(self) -> Iterator[Scale]:
         for rmin, rmax in self.as_array():
             yield Scale(rmin=rmin, rmax=rmax)
 
+    def modify(
+        self,
+        rmin: list[float] | float = DEFAULT.NotSet,
+        rmax: list[float] | float = DEFAULT.NotSet,
+        rweight: float | None = DEFAULT.NotSet,
+        rbin_num: int = DEFAULT.NotSet,
+    ) -> ScalesConfig:
+        return super().modify(rmin=rmin, rmax=rmax, rweight=rweight, rbin_num=rbin_num)
+
     def as_array(self) -> NDArray[np.float_]:
         """Obtain the scales cuts as array of shape (2, N)"""
         return np.atleast_2d(np.transpose([self.rmin, self.rmax]))
@@ -155,4 +164,4 @@ def dict_keys(self) -> list[str]:
             >>> [str(scale) for scale in ScalesConfig]
             ...
         """
-        return [str(scale) for scale in self]
+        return [str(scale) for scale in self]  # pragma: no cover
diff --git a/src/yaw/core/abc.py b/src/yaw/core/abc.py
index f3073aeb..50dbb264 100644
--- a/src/yaw/core/abc.py
+++ b/src/yaw/core/abc.py
@@ -4,7 +4,7 @@
 
 from __future__ import annotations
 
-from abc import ABC, abstractclassmethod, abstractmethod, abstractproperty
+from abc import ABC, abstractmethod
 from dataclasses import asdict
 from typing import TYPE_CHECKING, Any, Type, TypeVar
 
@@ -34,12 +34,14 @@
 class PatchedQuantity(ABC):
     """Base class for an object that has data organised in spatial patches."""
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def n_patches(self) -> int:
         """Get the number of spatial patches."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def patches(self) -> Indexer:
         """An :obj:`~yaw.core.containers.Indexer` attribute that supports
         iteration over the spatial patches or selecting a subset of the patches.
@@ -129,7 +131,8 @@ def closed(self) -> str:
         be: ``left``, ``right``, ``both``, ``neither``."""
         return self.get_binning().closed
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def bins(self) -> Indexer:
         """An :obj:`~yaw.core.containers.Indexer` attribute that supports
         iteration over the bins or selecting a subset of the bins.
@@ -227,7 +230,8 @@ def concatenate_bin_edges(*patched: BinnedQuantity) -> IntervalIndex:
 class HDFSerializable(ABC):
     """Base class for an object that can be serialised into a HDF5 file."""
 
-    @abstractclassmethod
+    @classmethod
+    @abstractmethod
     def from_hdf(cls: Type[_Thdf], source: h5py.Group) -> _Thdf:
         """Create a class instance by deserialising data from a HDF5 group.
 
diff --git a/src/yaw/core/containers.py b/src/yaw/core/containers.py
index a7bfa32a..26997496 100644
--- a/src/yaw/core/containers.py
+++ b/src/yaw/core/containers.py
@@ -174,19 +174,14 @@ def __str__(self) -> str:
         return f"{self.value:+.3g}+/-{self.error:.3g}"
 
     def __eq__(self, other: object) -> bool:
-        if isinstance(other, self.__class__):
-            if self.samples.shape != other.samples.shape:
-                return False
+        if isinstance(other, SampledValue):
             return (
-                self.method == other.method
+                self.samples.shape == other.samples.shape
+                and self.method == other.method
                 and self.value == other.value
                 and np.all(self.samples == other.samples)
             )
-        else:
-            return False
-
-    def __ne__(self, other: object) -> bool:
-        return not self.__eq__(other)
+        return NotImplemented
 
     @property
     def n_samples(self) -> int:
@@ -297,37 +292,36 @@ def __repr__(self) -> str:
 
     def __eq__(self, other: object) -> bool:
         if isinstance(other, self.__class__):
-            if self.samples.shape != other.samples.shape:
-                return False
             return (
-                self.method == other.method
+                self.samples.shape == other.samples.shape
+                and self.method == other.method
                 and np.all(self.data == other.data)
                 and np.all(self.samples == other.samples)
                 and (self.binning == other.binning).all()
             )
-        else:
-            return False
-
-    def __ne__(self, other: object) -> bool:
-        return not self.__eq__(other)
+        return NotImplemented
 
-    def __add__(self, other: _Tdata) -> _Tdata:
-        self.is_compatible(other, require=True)
-        return self.__class__(
-            binning=self.get_binning(),
-            data=self.data + other.data,
-            samples=self.samples + other.samples,
-            method=self.method,
-        )
+    def __add__(self, other: object) -> _Tdata:
+        if not isinstance(other, self.__class__):
+            self.is_compatible(other, require=True)
+            return self.__class__(
+                binning=self.get_binning(),
+                data=self.data + other.data,
+                samples=self.samples + other.samples,
+                method=self.method,
+            )
+        return NotImplemented
 
-    def __sub__(self, other: _Tdata) -> _Tdata:
-        self.is_compatible(other, require=True)
-        return self.__class__(
-            binning=self.get_binning(),
-            data=self.data - other.data,
-            samples=self.samples - other.samples,
-            method=self.method,
-        )
+    def __sub__(self, other: object) -> _Tdata:
+        if not isinstance(other, self.__class__):
+            self.is_compatible(other, require=True)
+            return self.__class__(
+                binning=self.get_binning(),
+                data=self.data - other.data,
+                samples=self.samples - other.samples,
+                method=self.method,
+            )
+        return NotImplemented
 
     @property
     def bins(self: _Tdata) -> Indexer[int | slice | Sequence, _Tdata]:
diff --git a/src/yaw/core/coordinates.py b/src/yaw/core/coordinates.py
index a6d114ce..d47d08bc 100644
--- a/src/yaw/core/coordinates.py
+++ b/src/yaw/core/coordinates.py
@@ -5,7 +5,7 @@
 
 from __future__ import annotations
 
-from abc import ABC, abstractclassmethod, abstractmethod, abstractproperty
+from abc import ABC, abstractmethod
 from collections.abc import Iterator, Sequence
 from functools import total_ordering
 from typing import TYPE_CHECKING, TypeVar
@@ -27,11 +27,13 @@ class Coordinate(ABC):
     def __init__(self, coords: dict[str, ArrayLike]) -> None:
         pass
 
-    @abstractclassmethod
+    @classmethod
+    @abstractmethod
     def from_array(cls, array) -> Coordinate:
         pass
 
-    @abstractclassmethod
+    @classmethod
+    @abstractmethod
     def from_coords(cls, coords: Sequence[Coordinate]) -> Coordinate:
         """Concatenate a sequence of coordinates into a new vector of
         coordates."""
@@ -40,7 +42,8 @@ def from_coords(cls, coords: Sequence[Coordinate]) -> Coordinate:
     def __repr__(self) -> str:
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def dim(self) -> tuple[str]:
         """A list of names of coordinates in the coordinate system."""
         pass
@@ -50,7 +53,8 @@ def ndim(self) -> int:
         """The number of coordinates in the coordinate system."""
         return len(self.dim)
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def values(self) -> NDArray[np.float_]:
         """The coordinate values cast into a numpy array with shape
         `(N, ndim)` or `(ndim,)` if there is only a single entry."""
@@ -160,8 +164,9 @@ def to_sky(self) -> CoordSky:
         r_d2 = np.sqrt(x * x + y * y)
         r_d3 = np.sqrt(x * x + y * y + z * z)
         # transform
-        ra = np.arccos(x / r_d2) * sgn(y) % (2.0 * np.pi)
-        ra[np.isnan(ra)] = 0.0
+        x_normed = np.ones_like(x)  # fallback for zero-division, arccos(1)=0.0
+        np.divide(x, r_d2, where=r_d2 > 0.0, out=x_normed)
+        ra = np.arccos(x_normed) * sgn(y) % (2.0 * np.pi)
         dec = np.arcsin(self.z / r_d3)
         return CoordSky(ra, dec)
 
@@ -334,18 +339,22 @@ def __repr__(self) -> str:
     def __format__(self, __format_spec: str) -> str:
         return self.values.__format__(__format_spec)
 
-    def __eq__(self, other: Distance) -> ArrayLike[np.bool_]:
-        return self.values == other.values
+    def __eq__(self, other: object) -> ArrayLike[np.bool_]:
+        if isinstance(other, self.__class__):
+            return self.values == other.values
+        return NotImplemented
 
     def __lt__(self, other: Distance) -> ArrayLike[np.bool_]:
-        return self.values < other.values
+        if isinstance(other, self.__class__):
+            return self.values < other.values
+        return NotImplemented
 
     @abstractmethod
-    def __add__(self: _Tdist, other: _Tdist) -> _Tdist:
+    def __add__(self, other: object) -> _Tdist:
         pass
 
     @abstractmethod
-    def __sub__(self: _Tdist, other: _Tdist) -> _Tdist:
+    def __sub__(self, other: object) -> _Tdist:
         pass
 
     def min(self) -> _Tdist:
@@ -372,13 +381,17 @@ def to_sky(self) -> DistSky:
 class Dist3D(Distance):
     """Implements a vector of Euclidean distances."""
 
-    def __add__(self, other: Dist3D) -> DistSky:
-        dist_sky = self.to_sky() + other.to_sky()
-        return dist_sky.to_3d()
+    def __add__(self, other: object) -> DistSky:
+        if isinstance(other, Dist3D):
+            dist_sky = self.to_sky() + other.to_sky()
+            return dist_sky.to_3d()
+        return NotImplemented
 
-    def __sub__(self, other: Dist3D) -> DistSky:
-        dist_sky = self.to_sky() - other.to_sky()
-        return dist_sky.to_3d()
+    def __sub__(self, other: object) -> DistSky:
+        if isinstance(other, Dist3D):
+            dist_sky = self.to_sky() - other.to_sky()
+            return dist_sky.to_3d()
+        return NotImplemented
 
     def to_3d(self) -> Dist3D:
         return self
@@ -392,11 +405,15 @@ def to_sky(self) -> DistSky:
 class DistSky(Distance):
     """Implements a vector of angular distances in radian."""
 
-    def __add__(self, other: DistSky) -> DistSky:
-        return DistSky(self.values + other.values)
+    def __add__(self, other: object) -> DistSky:
+        if isinstance(other, DistSky):
+            return DistSky(self.values + other.values)
+        return NotImplemented
 
-    def __sub__(self, other: DistSky) -> DistSky:
-        return DistSky(self.values - other.values)
+    def __sub__(self, other: object) -> DistSky:
+        if isinstance(other, DistSky):
+            return DistSky(self.values - other.values)
+        return NotImplemented
 
     def to_3d(self) -> Dist3D:
         return Dist3D(2.0 * np.sin(self.values / 2.0))
diff --git a/src/yaw/core/parallel.py b/src/yaw/core/parallel.py
index 19c48d4e..110649bf 100644
--- a/src/yaw/core/parallel.py
+++ b/src/yaw/core/parallel.py
@@ -279,7 +279,7 @@ def _init_pool(
                 all_initargs.append(None)
             else:
                 raise TypeError("'initargs' must be an iterable")
-        logger.debug(f"running {self.n_jobs()} jobs on {self.n_threads()} threads")
+        logger.debug("running %i jobs on %i threads", self.n_jobs(), self.n_threads())
         return multiprocessing.Pool(
             initializer=_threadinit, initargs=all_initargs, processes=self._num_threads
         )
diff --git a/src/yaw/correlation/corrfuncs.py b/src/yaw/correlation/corrfuncs.py
index f537ba7b..b989af01 100644
--- a/src/yaw/correlation/corrfuncs.py
+++ b/src/yaw/correlation/corrfuncs.py
@@ -128,7 +128,7 @@ def from_files(cls: Type[_Tdata], path_prefix: TypePathStr) -> _Tdata:
             :obj:`CorrData`
         """
         name = cls.__name__.lower()[:-4]
-        logger.debug(f"reading {name} data from '{path_prefix}.*'")
+        logger.debug("reading %s data from '%s.*'", name, path_prefix)
         # load data and errors
         ext = "dat"
         data_error = np.loadtxt(f"{path_prefix}.{ext}")
@@ -212,7 +212,7 @@ def to_files(self, path_prefix: TypePathStr) -> None:
                 The base name of the output files without any file extension.
         """
         name = self.__class__.__name__.lower()[:-4]
-        logger.info(f"writing {name} data to '{path_prefix}.*'")
+        logger.info("writing %s data to '%s.*'", name, path_prefix)
         PREC = 10
         DELIM = " "
 
@@ -306,7 +306,7 @@ def plot(
         plot_kwargs: dict[str, Any] | None = None,
         zero_line: bool = False,
         scale_by_dz: bool = False,
-    ) -> Axis:  # pragma: no cover
+    ) -> Axis:
         """Create a plot of the correlation data as a function of redshift.
 
         Create a new axis or plot to an existing one, add x-axis offsets, if
@@ -337,8 +337,8 @@ def plot(
                 :obj:`dz`.
         """
         x = self.mids + xoffset
-        y = self.data
-        yerr = self.error
+        y = self.data.astype(np.float_)
+        yerr = self.error.astype(np.float_)
         if scale_by_dz:
             y *= self.dz
             yerr *= self.dz
@@ -493,7 +493,8 @@ def __post_init__(self) -> None:
                 continue
             try:
                 self.dd.is_compatible(pairs, require=True)
-            except ValueError as e:
+                assert self.dd.n_patches == pairs.n_patches
+            except (ValueError, AssertionError) as e:
                 raise ValueError(
                     f"pair counts '{kind}' and 'dd' are not compatible"
                 ) from e
@@ -505,48 +506,52 @@ def __repr__(self) -> str:
         other = f"n_patches={self.n_patches}"
         return f"{string}, {pairs}, {other})"
 
-    def __eq__(self, other) -> bool:
-        if not isinstance(other, self.__class__):
-            return False
-        for cfield in fields(self):
-            kind = cfield.name
-            if getattr(self, kind) != getattr(other, kind):
-                return False
-        return True
-
-    def __neq__(self, other) -> bool:
-        return not self == other
-
-    def __add__(self, other: CorrFunc) -> CorrFunc:
-        # check that the pair counts are set consistently
-        kinds = []
-        for cfield in fields(self):
-            kind = cfield.name
-            self_set = getattr(self, kind) is not None
-            other_set = getattr(other, kind) is not None
-            if (self_set and not other_set) or (not self_set and other_set):
-                raise ValueError(f"pair counts for '{kind}' not set for both operands")
-            elif self_set and other_set:
-                kinds.append(kind)
-
-        kwargs = {kind: getattr(self, kind) + getattr(other, kind) for kind in kinds}
-        return self.__class__(**kwargs)
-
-    def __radd__(self, other: CorrFunc | int | float) -> CorrFunc:
-        if other == 0:
+    def __eq__(self, other: object) -> bool:
+        if isinstance(other, self.__class__):
+            for cfield in fields(self):
+                kind = cfield.name
+                if getattr(self, kind) != getattr(other, kind):
+                    return False
+            return True
+        return NotImplemented
+
+    def __add__(self, other: object) -> CorrFunc:
+        if isinstance(other, self.__class__):
+            # check that the pair counts are set consistently
+            kinds = []
+            for cfield in fields(self):
+                kind = cfield.name
+                self_set = getattr(self, kind) is not None
+                other_set = getattr(other, kind) is not None
+                if (self_set and not other_set) or (not self_set and other_set):
+                    raise ValueError(
+                        f"pair counts for '{kind}' not set for both operands"
+                    )
+                elif self_set and other_set:
+                    kinds.append(kind)
+
+            kwargs = {
+                kind: getattr(self, kind) + getattr(other, kind) for kind in kinds
+            }
+            return self.__class__(**kwargs)
+        return NotImplemented
+
+    def __radd__(self, other: object) -> CorrFunc:
+        if np.isscalar(other) and other == 0:
             return self
-        else:
-            return self.__add__(other)
-
-    def __mul__(self, other: np.number) -> CorrFunc:
-        # check that the pair counts are set consistently
-        kwargs = {}
-        for cfield in fields(self):
-            kind = cfield.name
-            counts = getattr(self, kind)
-            if counts is not None:
-                kwargs[kind] = counts * other
-        return self.__class__(**kwargs)
+        return other.__add__(self)
+
+    def __mul__(self, other: object) -> CorrFunc:
+        if np.isscalar(other) and not isinstance(other, (bool, np.bool_)):
+            # check that the pair counts are set consistently
+            kwargs = {}
+            for cfield in fields(self):
+                kind = cfield.name
+                counts = getattr(self, kind)
+                if counts is not None:
+                    kwargs[kind] = counts * other
+            return self.__class__(**kwargs)
+        return NotImplemented
 
     @property
     def auto(self) -> bool:
@@ -589,7 +594,7 @@ def get_binning(self) -> IntervalIndex:
     def n_patches(self) -> int:
         return self.dd.n_patches
 
-    def is_compatible(self, other: CorrFunc, require: bool = True) -> bool:
+    def is_compatible(self, other: CorrFunc, require: bool = False) -> bool:
         """Check whether this instance is compatible with another instance.
 
         Ensures that the redshift binning and the number of patches are
@@ -623,8 +628,6 @@ def estimators(self) -> dict[str, CorrelationEstimator]:
         available = set()
         # iterate all dataclass attributes that are in __init__
         for attr in fields(self):
-            if not attr.init:
-                continue
             if getattr(self, attr.name) is not None:
                 available.add(cts_from_code(attr.name))
         # check which estimators are supported
@@ -635,7 +638,7 @@ def estimators(self) -> dict[str, CorrelationEstimator]:
         return estimators
 
     def _check_and_select_estimator(
-        self, estimator: str | None
+        self, estimator: str | None = None
     ) -> type[CorrelationEstimator]:
         options = self.estimators
         if estimator is None:
@@ -655,17 +658,13 @@ def _check_and_select_estimator(
             # determine which pair counts are missing
             for attr in fields(self):
                 name = attr.name
-                if not attr.init:
-                    continue
                 cts = cts_from_code(name)
                 if getattr(self, name) is None and cts in est_class.requires:
                     raise EstimatorError(f"estimator requires {name}")
-            else:
-                raise RuntimeError()
         # select the correct estimator
         cls = options[estimator]
         logger.debug(
-            f"selecting estimator '{cls.short}' from {'/'.join(self.estimators)}"
+            "selecting estimator '%s' from %s", cls.short, "/".join(self.estimators)
         )
         return cls
 
@@ -685,7 +684,7 @@ def get(self, *args, **kwargs):
         .. deprecated:: 2.3.1
             Renamed to :meth:`sample`.
         """
-        return self.sample(*args, **kwargs)
+        return self.sample(*args, **kwargs)  # pragma: no cover
 
     def sample(
         self,
@@ -718,7 +717,7 @@ def sample(
         if config is None:
             config = ResamplingConfig()
         est_fun = self._check_and_select_estimator(estimator)
-        logger.debug(f"computing correlation and {config.method} samples")
+        logger.debug("computing correlation and %s samples", config.method)
         # get the pair counts for the required terms (DD, maybe DR and/or RR)
         required_data = {}
         required_samples = {}
@@ -779,12 +778,12 @@ def to_hdf(self, dest: h5py.File | h5py.Group) -> None:
 
     @classmethod
     def from_file(cls, path: TypePathStr) -> CorrFunc:
-        logger.debug(f"reading pair counts from '{path}'")
+        logger.debug("reading pair counts from '%s'", path)
         with h5py.File(str(path)) as f:
             return cls.from_hdf(f)
 
     def to_file(self, path: TypePathStr) -> None:
-        logger.info(f"writing pair counts to '{path}'")
+        logger.info("writing pair counts to '%s'", path)
         with h5py.File(str(path), mode="w") as f:
             self.to_hdf(f)
 
@@ -928,8 +927,10 @@ def autocorrelate(
     _check_patch_centers([data, random])
     scales = config.scales.as_array()
     logger.info(
-        f"running autocorrelation ({len(scales)} scales, "
-        f"{scales.min():.0f}<r<={scales.max():.0f}kpc)"
+        "running autocorrelation (%i scales, %.0f<r<=%.0fkpc)",
+        len(scales),
+        scales.min(),
+        scales.max(),
     )
     if linkage is None:
         linkage = PatchLinkage.from_setup(config, random)
@@ -1023,8 +1024,10 @@ def crosscorrelate(
 
     scales = config.scales.as_array()
     logger.info(
-        f"running crosscorrelation ({len(scales)} scales, "
-        f"{scales.min():.0f}<r<={scales.max():.0f}kpc)"
+        "running crosscorrelation (%i scales, %.0f<r<=%.0fkpc)",
+        len(scales),
+        scales.min(),
+        scales.max(),
     )
     if linkage is None:
         linkage = PatchLinkage.from_setup(config, unknown)
diff --git a/src/yaw/correlation/estimators.py b/src/yaw/correlation/estimators.py
index 243aa8e2..8ce14f40 100644
--- a/src/yaw/correlation/estimators.py
+++ b/src/yaw/correlation/estimators.py
@@ -23,7 +23,7 @@
 from __future__ import annotations
 
 import warnings
-from abc import ABC, abstractclassmethod, abstractproperty
+from abc import ABC, abstractmethod
 from typing import TYPE_CHECKING
 
 import numpy as np
@@ -47,7 +47,8 @@ class EstimatorError(Exception):
 class Cts(ABC):
     """Base class to symbolically represent pair counts."""
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def _hash(self) -> int:
         """Used for comparison.
 
@@ -55,7 +56,8 @@ def _hash(self) -> int:
         :obj:`CtsMix`."""
         pass
 
-    @abstractproperty
+    @property
+    @abstractmethod
     def _str(self) -> str:
         pass
 
@@ -69,11 +71,11 @@ def __hash__(self) -> int:
         return self._hash
 
     def __eq__(self, other: object) -> bool:
-        if not isinstance(other, Cts):
-            return False
-        var1 = set(self._str.split("_"))
-        var2 = set(other._str.split("_"))
-        return not var1.isdisjoint(var2)
+        if isinstance(other, Cts):
+            var1 = set(self._str.split("_"))
+            var2 = set(other._str.split("_"))
+            return not var1.isdisjoint(var2)
+        return NotImplemented
 
 
 class CtsDD(Cts):
@@ -158,7 +160,8 @@ def _warn_enum_zero(cls, counts: NDArray):
         if np.any(counts == 0.0):
             warnings.warn(f"estimator {cls.short} encontered zeros in enumerator")
 
-    @abstractclassmethod
+    @classmethod
+    @abstractmethod
     def eval(
         cls, *, dd: NDArray, dr: NDArray, rr: NDArray, rd: NDArray | None = None
     ) -> NDArray:
diff --git a/src/yaw/correlation/paircounts.py b/src/yaw/correlation/paircounts.py
index 1360a313..52d9c4a0 100644
--- a/src/yaw/correlation/paircounts.py
+++ b/src/yaw/correlation/paircounts.py
@@ -17,7 +17,7 @@
 from collections.abc import Iterable, Sequence
 from dataclasses import dataclass
 from itertools import accumulate
-from typing import TYPE_CHECKING, NoReturn, Union
+from typing import TYPE_CHECKING, NoReturn, Type, Union
 
 try:  # pragma: no cover
     from typing import TypeAlias
@@ -53,6 +53,12 @@
 TypeIndex: TypeAlias = Union[int, slice, Sequence]
 
 
+def sequence_require_type(items: Sequence, class_or_inst: Type | object) -> None:
+    for item in items:
+        if not isinstance(item, class_or_inst):
+            raise TypeError(f"invalid type '{type(item)}' for concatenation")
+
+
 def check_mergable(patched_arrays: Sequence[PatchedArray], *, patches: bool) -> None:
     """Check if two instaces of PatchedArray can be merged along the patch
     or binning axis.
@@ -71,11 +77,8 @@ def check_mergable(patched_arrays: Sequence[PatchedArray], *, patches: bool) ->
             raise ValueError("cannot merge mixed cross- and autocorrelations")
         if patches:
             reference.is_compatible(patched, require=True)
-        else:
-            if reference.auto != patched.auto:
-                raise ValueError("cannot merge mixed cross- and autocorrelations")
-            if reference.n_patches != patched.n_patches:
-                raise ValueError("cannot merge, patch numbers do not match")
+        elif reference.n_patches != patched.n_patches:
+            raise ValueError("cannot merge, patch numbers do not match")
 
 
 def binning_from_hdf(source: h5py.Group) -> IntervalIndex:
@@ -122,18 +125,8 @@ def __repr__(self) -> str:
         return f"{string}, {shape=})"
 
     @abstractmethod
-    def __eq__(self, other) -> bool:
-        if isinstance(other, self.__class__):
-            if self.n_bins != other.n_bins:
-                return False
-            elif self.n_patches != other.n_patches:
-                return False
-            elif not (self.get_binning() == other.get_binning()).all():
-                return False
-        return True
-
-    def __neq__(self, other) -> bool:
-        return not self == other
+    def __eq__(self, other: object) -> bool:
+        pass
 
     @property
     def dtype(self) -> DTypeLike:
@@ -193,7 +186,7 @@ def get_sum(self, *args, **kwargs):
         .. deprecated:: 2.3.1
             Renamed to :meth:`sample_sum`.
         """
-        return self.sample_sum(*args, **kwargs)
+        return self.sample_sum(*args, **kwargs)  # pragma: no cover
 
     def sample_sum(self, config: ResamplingConfig | None = None) -> SampledData:
         """Compute the sum of the data over all patches and samples thereof.
@@ -328,14 +321,17 @@ def __init__(
         self.totals2 = totals2
         self.auto = auto
 
-    def __eq__(self, other) -> bool:
-        if not super().__eq__(other):
-            return False  # checks type
-        return (
-            np.all(self.totals1 == other.totals1)
-            and np.all(self.totals2 == other.totals2)
-            and self.auto == other.auto
-        )
+    def __eq__(self, other: object) -> bool:
+        if isinstance(other, self.__class__):
+            return (
+                self.n_bins == other.n_bins
+                and self.n_patches == other.n_patches
+                and (self.get_binning() == other.get_binning()).all()
+                and np.all(self.totals1 == other.totals1)
+                and np.all(self.totals2 == other.totals2)
+                and self.auto == other.auto
+            )
+        return NotImplemented
 
     def as_array(self) -> NDArray:
         return np.einsum("i...,j...->ij...", self.totals1, self.totals2)
@@ -394,6 +390,7 @@ def to_hdf(self, dest: h5py.Group) -> None:
         dest.create_dataset("auto", data=self.auto)
 
     def concatenate_patches(self, *data: PatchedTotal) -> PatchedTotal:
+        sequence_require_type(data, self.__class__)
         all_totals: list[PatchedTotal] = [self, *data]
         check_mergable(all_totals, patches=True)
         return self.__class__(
@@ -404,6 +401,7 @@ def concatenate_patches(self, *data: PatchedTotal) -> PatchedTotal:
         )
 
     def concatenate_bins(self, *data: PatchedTotal) -> PatchedTotal:
+        sequence_require_type(data, self.__class__)
         all_totals: list[PatchedTotal] = [self, *data]
         check_mergable(all_totals, patches=False)
         binning = concatenate_bin_edges(*all_totals)
@@ -640,27 +638,38 @@ def zeros(
         counts = np.zeros((n_patches, n_patches, len(binning)), dtype=dtype)
         return cls(binning, counts, auto=auto)
 
-    def __eq__(self, other) -> bool:
-        if not super().__eq__(other):
-            return False  # checks type
-        return np.all(self.counts == other.counts) and (self.auto == other.auto)
+    def __eq__(self, other: object) -> bool:
+        if isinstance(other, self.__class__):
+            return (
+                self.n_bins == other.n_bins
+                and self.n_patches == other.n_patches
+                and (self.get_binning() == other.get_binning()).all()
+                and np.all(self.counts == other.counts)
+                and (self.auto == other.auto)
+            )
+        return NotImplemented
 
-    def __add__(self, other: PatchedCount) -> PatchedCount:
-        self.is_compatible(other, require=True)
-        if self.n_patches != other.n_patches:
-            raise ValueError("number of patches does not agree")
-        return self.__class__(
-            self.get_binning(), self.counts + other.counts, auto=self.auto
-        )
+    def __add__(self, other: object) -> PatchedCount:
+        if isinstance(other, self.__class__):
+            self.is_compatible(other, require=True)
+            if self.n_patches != other.n_patches:
+                raise ValueError("number of patches does not agree")
+            return self.__class__(
+                self.get_binning(), self.counts + other.counts, auto=self.auto
+            )
+        return NotImplemented
 
-    def __radd__(self, other: PatchedCount | int | float) -> PatchedCount:
-        if other == 0:
+    def __radd__(self, other: object) -> PatchedCount:
+        if np.isscalar(other) and other == 0:
             return self
-        else:
-            return self.__add__(other)
+        return other.__add__(self)
 
-    def __mul__(self, other: np.number) -> PatchedCount:
-        return self.__class__(self.get_binning(), self.counts * other, auto=self.auto)
+    def __mul__(self, other: object) -> PatchedCount:
+        if np.isscalar(other) and not isinstance(other, (bool, np.bool_)):
+            return self.__class__(
+                self.get_binning(), self.counts * other, auto=self.auto
+            )
+        return NotImplemented
 
     def set_measurement(self, key: PatchIDs | tuple[int, int], item: NDArray):
         """Set the counts value in all redshift bins for a pair of patch
@@ -783,22 +792,27 @@ def to_hdf(self, dest: h5py.Group) -> None:
         dest.create_dataset("auto", data=self.auto)
 
     def concatenate_patches(self, *data: PatchedCount) -> PatchedCount:
+        sequence_require_type(data, self.__class__)
         all_counts: list[PatchedCount] = [self, *data]
         check_mergable(all_counts, patches=True)
         offsets = patch_idx_offset(all_counts)
+        n_patches = [count.n_patches for count in all_counts]
         merged = self.__class__.zeros(
             binning=self.get_binning(),
-            n_patches=sum(count.n_patches for count in all_counts),
+            n_patches=sum(n_patches),
             auto=self.auto,
         )
         # insert the blocks of counts into the merged counts array
         loc = 0
-        for count, offset in zip(all_counts, offsets):
-            merged.counts[loc : loc + offset, loc : loc + offset] = count.counts
+        for count, offset, n in zip(all_counts, offsets, n_patches):
+            i_start = offset
+            i_end = i_start + n
+            merged.counts[i_start:i_end, i_start:i_end] = count.counts
             loc += offset
         return merged
 
     def concatenate_bins(self, *data: PatchedCount) -> PatchedCount:
+        sequence_require_type(data, self.__class__)
         all_counts: list[PatchedCount] = [self, *data]
         check_mergable(all_counts, patches=False)
         binning = concatenate_bin_edges(*all_counts)
@@ -952,33 +966,30 @@ def __repr__(self) -> str:
         n_patches = self.n_patches
         return f"{string}, {n_patches=})"
 
-    def __eq__(self, other) -> bool:
+    def __eq__(self, other: object) -> bool:
         if isinstance(other, self.__class__):
-            return np.all(self.count == other.count) and np.all(
-                self.total == other.total
-            )
-        else:
-            return False
-
-    def __neq__(self, other) -> bool:
-        return not self == other
-
-    def __add__(self, other: NormalisedCounts) -> NormalisedCounts:
-        count = self.count + other.count
-        if np.any(self.total.totals1 != other.total.totals1) or np.any(
-            self.total.totals2 != other.total.totals2
-        ):
-            raise ValueError("total number of objects do not agree for operands")
-        return self.__class__(count, self.total)
+            return self.count == other.count and self.total == other.total
+        return NotImplemented
 
-    def __radd__(self, other: NormalisedCounts | int | float) -> NormalisedCounts:
-        if other == 0:
+    def __add__(self, other: object) -> NormalisedCounts:
+        if isinstance(other, self.__class__):
+            count = self.count + other.count
+            if np.any(self.total.totals1 != other.total.totals1) or np.any(
+                self.total.totals2 != other.total.totals2
+            ):
+                raise ValueError("total number of objects do not agree for operands")
+            return self.__class__(count, self.total)
+        return NotImplemented
+
+    def __radd__(self, other: object) -> NormalisedCounts:
+        if np.isscalar(other) and other == 0:
             return self
-        else:
-            return self.__add__(other)
+        return other.__add__(self)
 
-    def __mul__(self, other: np.number) -> NormalisedCounts:
-        return self.__class__(self.count * other, self.total)
+    def __mul__(self, other: object) -> NormalisedCounts:
+        if np.isscalar(other) and not isinstance(other, (bool, np.bool_)):
+            return self.__class__(self.count * other, self.total)
+        return NotImplemented
 
     @property
     def auto(self) -> bool:
@@ -1054,6 +1065,7 @@ def to_hdf(self, dest: h5py.Group) -> None:
         self.total.to_hdf(group)
 
     def concatenate_patches(self, *pcounts: NormalisedCounts) -> NormalisedCounts:
+        sequence_require_type(pcounts, self.__class__)
         counts = [pc.count for pc in pcounts]
         totals = [pc.total for pc in pcounts]
         return self.__class__(
@@ -1062,6 +1074,7 @@ def concatenate_patches(self, *pcounts: NormalisedCounts) -> NormalisedCounts:
         )
 
     def concatenate_bins(self, *pcounts: NormalisedCounts) -> NormalisedCounts:
+        sequence_require_type(pcounts, self.__class__)
         counts = [pc.count for pc in pcounts]
         totals = [pc.total for pc in pcounts]
         return self.__class__(
diff --git a/src/yaw/deprecated.py b/src/yaw/deprecated.py
deleted file mode 100644
index a7932af7..00000000
--- a/src/yaw/deprecated.py
+++ /dev/null
@@ -1,53 +0,0 @@
-from __future__ import annotations
-
-from deprecated import deprecated
-
-from yaw.correlation import CorrData, CorrFunc
-from yaw.correlation.paircounts import NormalisedCounts
-from yaw.redshifts import HistData
-
-__all__ = ["CorrelationData", "CorrelationFunction", "HistogramData", "PairCountResult"]
-
-
-@deprecated(reason="renamed to yaw.CorrData", action="module", version="2.3.2")
-class CorrelationData(CorrData):
-    """
-    .. deprecated:: 2.3.2
-        Renamed to :meth:`yaw.CorrData`.
-    """
-
-    pass
-
-
-@deprecated(reason="renamed to yaw.CorrFunc", action="module", version="2.3.2")
-class CorrelationFunction(CorrFunc):
-    """
-    .. deprecated:: 2.3.2
-        Renamed to :meth:`yaw.CorrFunc`.
-    """
-
-    pass
-
-
-@deprecated(reason="renamed to yaw.HistData", action="module", version="2.3.2")
-class HistogramData(HistData):
-    """
-    .. deprecated:: 2.3.2
-        Renamed to :meth:`yaw.HistData`.
-    """
-
-    pass
-
-
-@deprecated(
-    reason="renamed to yaw.correlation.paircounts.NormalisedCounts",
-    action="module",
-    version="2.3.2",
-)
-class PairCountResult(NormalisedCounts):
-    """
-    .. deprecated:: 2.3.2
-        Renamed to :obj:`yaw.correlation.paircounts.NormalisedCounts`.
-    """
-
-    pass
diff --git a/src/yaw/deprecated/__init__.py b/src/yaw/deprecated/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/src/yaw/deprecated/config/__init__.py b/src/yaw/deprecated/config/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/src/yaw/deprecated/config/binning.py b/src/yaw/deprecated/config/binning.py
new file mode 100644
index 00000000..3104093f
--- /dev/null
+++ b/src/yaw/deprecated/config/binning.py
@@ -0,0 +1,266 @@
+from __future__ import annotations
+
+import warnings
+from dataclasses import dataclass, field
+from typing import TYPE_CHECKING, Any
+
+import numpy as np
+from deprecated import deprecated
+
+from yaw.config import OPTIONS
+from yaw.config import default as DEFAULT
+from yaw.config.utils import ConfigError
+from yaw.core.abc import DictRepresentation
+from yaw.core.cosmology import BinFactory, TypeCosmology
+from yaw.core.docs import Parameter
+from yaw.core.math import array_equal
+
+if TYPE_CHECKING:  # pragma: no cover
+    from numpy.typing import NDArray
+
+__all__ = ["AutoBinningConfig", "ManualBinningConfig", "make_binning_config"]
+
+
+class BaseBinningConfig(DictRepresentation):
+    """Base class for redshift binning configuration."""
+
+    zbins: NDArray[np.float_]
+    """Edges of redshift bins."""
+    method: str
+    """Method used to create redshift binning, ``manual`` or either of
+    :obj:`~yaw.config.options.Options.binning`."""
+
+    def __repr__(self) -> str:
+        name = self.__class__.__name__
+        zbin_num = self.zbin_num
+        z = f"{self.zmin:.3f}...{self.zmax:.3f}"
+        method = self.method
+        return f"{name}({zbin_num=}, {z=}, {method=})"
+
+
+@deprecated(
+    reason="merged into yaw.config.BinningConfig", action="module", version="2.5.5"
+)
+@dataclass(frozen=True, repr=False)
+class ManualBinningConfig(BaseBinningConfig):
+    """Configuration that specifies a manual redshift binning.
+
+    Args:
+        zbins (:obj:`NDArray`):
+            Edges of redshift bins, must increase monotonically.
+    """
+
+    zbins: NDArray[np.float_] = field(
+        metadata=Parameter(
+            type=float,
+            nargs="*",
+            help="list of custom redshift bin edges, if method is set to 'manual'",
+        )
+    )
+
+    def __post_init__(self) -> None:
+        if len(self.zbins) < 2:
+            raise ConfigError("'zbins' must have at least two edges")
+        BinFactory.check(self.zbins)
+        object.__setattr__(self, "zbins", np.asarray(self.zbins))
+
+    def __eq__(self, other: ManualBinningConfig) -> bool:
+        if not array_equal(self.zbins, other.zbins):
+            return False
+        return True
+
+    @property
+    def method(self) -> str:
+        """Method used to create redshift binning, always ``manual``."""
+        return "manual"
+
+    @property
+    def zmin(self) -> float:
+        """Lowest redshift bin edge."""
+        return float(self.zbins[0])
+
+    @property
+    def zmax(self) -> float:
+        """Highest redshift bin edge."""
+        return float(self.zbins[-1])
+
+    @property
+    def zbin_num(self) -> int:
+        """Number of redshift bins."""
+        return len(self.zbins) - 1
+
+    @classmethod
+    def from_dict(cls, the_dict: dict[str, Any], **kwargs) -> ManualBinningConfig:
+        new_dict = {key: val for key, val in the_dict.items()}
+        zbins = new_dict.pop("zbins")
+        return cls(np.asarray(zbins), **new_dict)
+
+    def to_dict(self) -> dict[str, Any]:
+        return dict(method=self.method, zbins=self.zbins.tolist())
+
+
+@deprecated(
+    reason="merged into yaw.config.BinningConfig", action="module", version="2.5.5"
+)
+@dataclass(frozen=True, repr=False)
+class AutoBinningConfig(BaseBinningConfig):
+    """Configuration that generates a redshift binning.
+
+    To generate a redshift binning use the :meth:`generate` method.
+
+    Args:
+        zbins (:obj:`NDArray`):
+            Edges of redshift bins, must increase monotonically.
+        method (:obj:`str`):
+            Method used to create redshift binning, either of
+            :obj:`~yaw.config.options.Options.binning`.
+    """
+
+    zbins: NDArray[np.float_]
+    method: str = field(
+        default=DEFAULT.Configuration.binning.method,
+        metadata=Parameter(
+            choices=OPTIONS.binning,
+            help="redshift binning method, 'logspace' means equal size in log(1+z)",
+            default_text="(default: %(default)s)",
+        ),
+    )
+    """Method used to create redshift binning, see
+    :obj:`~yaw.config.options.Options.binning`."""
+    zmin: float = field(
+        init=False,
+        metadata=Parameter(
+            type=float,
+            help="lower redshift limit",
+            default_text="(default: %(default)s)",
+        ),
+    )
+    """Lowest redshift bin edge."""
+    zmax: float = field(
+        init=False,
+        metadata=Parameter(
+            type=float,
+            help="upper redshift limit",
+            default_text="(default: %(default)s)",
+        ),
+    )
+    """Highest redshift bin edge."""
+    zbin_num: int = field(
+        default=DEFAULT.Binning.zbin_num,
+        init=False,
+        metadata=Parameter(
+            type=int,
+            help="number of redshift bins",
+            default_text="(default: %(default)s)",
+        ),
+    )
+    """Number of redshift bins."""
+
+    def __post_init__(self) -> None:
+        object.__setattr__(self, "zmin", float(self.zbins[0]))
+        object.__setattr__(self, "zmax", float(self.zbins[-1]))
+        object.__setattr__(self, "zbin_num", len(self.zbins) - 1)
+
+    @classmethod
+    def generate(
+        cls,
+        zmin: float,
+        zmax: float,
+        zbin_num: int = DEFAULT.Binning.zbin_num,
+        method: str = DEFAULT.Binning.method,
+        cosmology: TypeCosmology | None = None,
+    ) -> AutoBinningConfig:
+        """Generate a new redshift binning configuration.
+
+        Generates a specified number of bins between a lower and upper redshift
+        limit. The spacing of the bins depends the generation method, the
+        default is a linear spacing.
+
+        Args:
+            zmin (:obj:`float`):
+                Minimum redshift, lowest redshift edge.
+            zmax (:obj:`float`):
+                Maximum redshift, highest redshift edge.
+            zbin_num (:obj:`int`, optional):
+                Number of redshift bins to generate.
+            method (:obj:`str`, optional):
+                Method used to create redshift binning, for a list of valid
+                options and their description see
+                :obj:`~yaw.config.options.Options.binning`.
+            cosmology (:obj:`astropy.cosmology.FLRW`, :obj:`~yaw.core.cosmology.CustomCosmology`, optional):
+                Cosmological model used for distance calculations. For a custom
+                model, refer to :mod:`yaw.core.cosmology`.
+
+        Returns:
+            :obj:`AutoBinningConfig`
+        """
+        if not isinstance(method, str):
+            raise ValueError("'method' must a string")
+        zbins = BinFactory(zmin, zmax, zbin_num, cosmology).get(method)
+        return cls(zbins, method)
+
+    def __eq__(self, other: ManualBinningConfig) -> bool:
+        if self.method != other.method:
+            return False
+        if not array_equal(self.zbins, other.zbins):
+            return False
+        return True
+
+    @classmethod
+    def from_dict(
+        cls, the_dict: dict[str, Any], cosmology: TypeCosmology | None = None
+    ) -> AutoBinningConfig:
+        """Create a class instance from a dictionary representation of the
+        minimally required data.
+
+        Args:
+            the_dict (:obj:`dict`):
+                Dictionary containing the data.
+            cosmology (:obj:`astropy.cosmology.FLRW`, :obj:`~yaw.core.cosmology.CustomCosmology`, optional):
+                Cosmological model used for distance calculations. For a custom
+                model, refer to :mod:`yaw.core.cosmology`.
+        """
+        return cls.generate(**the_dict, cosmology=cosmology)
+
+    def to_dict(self) -> dict[str, Any]:
+        return dict(
+            zmin=self.zmin, zmax=self.zmax, zbin_num=self.zbin_num, method=self.method
+        )
+
+
+@deprecated(
+    reason="use yaw.config.BinningConfig.create() instead",
+    action="module",
+    version="2.5.5",
+)
+def make_binning_config(
+    *,
+    cosmology: TypeCosmology | str | None = None,
+    zmin: float | None = None,
+    zmax: float | None = None,
+    zbin_num: int = DEFAULT.Binning.zbin_num,
+    method: str = DEFAULT.Binning.method,
+    zbins: NDArray[np.float_] | None = None,
+) -> ManualBinningConfig | AutoBinningConfig:
+    """
+    Helper function to construct a binning configuration.
+
+    The ``cosmology`` argument is always required. If redshift bins (``zbins``)
+    are provided, a :obj:`ManualBinningConfig` is returned, otherwise an
+    :obj:`AutoBinningConfig`. Issues a warning if any argument is set but
+    ignored by the returned configuration instance.
+    """
+    auto_args_set = (zmin is not None, zmax is not None)
+    manual_args_set = (zbins is not None,)
+    if not all(manual_args_set) and not all(auto_args_set):
+        raise ConfigError("either 'zbins' or 'zmin' and 'zmax' are required")
+    elif all(auto_args_set):
+        if all(manual_args_set):
+            warnings.warn(
+                "'zbins' set but ignored since 'zmin' and 'zmax' are provided"
+            )
+        return AutoBinningConfig.generate(
+            zmin=zmin, zmax=zmax, zbin_num=zbin_num, method=method, cosmology=cosmology
+        )
+    else:
+        return ManualBinningConfig(zbins)
diff --git a/src/yaw/deprecated/correlation/__init__.py b/src/yaw/deprecated/correlation/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/src/yaw/deprecated/correlation/corrfuncs.py b/src/yaw/deprecated/correlation/corrfuncs.py
new file mode 100644
index 00000000..c78d9803
--- /dev/null
+++ b/src/yaw/deprecated/correlation/corrfuncs.py
@@ -0,0 +1,27 @@
+from __future__ import annotations
+
+from deprecated import deprecated
+
+from yaw.correlation import CorrData, CorrFunc
+
+__all__ = ["CorrelationData", "CorrelationFunction"]
+
+
+@deprecated(reason="renamed to yaw.CorrData", action="module", version="2.3.2")
+class CorrelationData(CorrData):
+    """
+    .. deprecated:: 2.3.2
+        Renamed to :meth:`yaw.CorrData`.
+    """
+
+    pass
+
+
+@deprecated(reason="renamed to yaw.CorrFunc", action="module", version="2.3.2")
+class CorrelationFunction(CorrFunc):
+    """
+    .. deprecated:: 2.3.2
+        Renamed to :meth:`yaw.CorrFunc`.
+    """
+
+    pass
diff --git a/src/yaw/deprecated/correlation/paircounts.py b/src/yaw/deprecated/correlation/paircounts.py
new file mode 100644
index 00000000..5e34a097
--- /dev/null
+++ b/src/yaw/deprecated/correlation/paircounts.py
@@ -0,0 +1,21 @@
+from __future__ import annotations
+
+from deprecated import deprecated
+
+from yaw.correlation.paircounts import NormalisedCounts
+
+__all__ = ["PairCountResult"]
+
+
+@deprecated(
+    reason="renamed to yaw.correlation.paircounts.NormalisedCounts",
+    action="module",
+    version="2.3.2",
+)
+class PairCountResult(NormalisedCounts):
+    """
+    .. deprecated:: 2.3.2
+        Renamed to :obj:`yaw.correlation.paircounts.NormalisedCounts`.
+    """
+
+    pass
diff --git a/src/yaw/deprecated/redshifts.py b/src/yaw/deprecated/redshifts.py
new file mode 100644
index 00000000..ee7178a6
--- /dev/null
+++ b/src/yaw/deprecated/redshifts.py
@@ -0,0 +1,17 @@
+from __future__ import annotations
+
+from deprecated import deprecated
+
+from yaw.redshifts import HistData
+
+__all__ = ["HistogramData"]
+
+
+@deprecated(reason="renamed to yaw.HistData", action="module", version="2.3.2")
+class HistogramData(HistData):
+    """
+    .. deprecated:: 2.3.2
+        Renamed to :meth:`yaw.HistData`.
+    """
+
+    pass
diff --git a/src/yaw/redshifts.py b/src/yaw/redshifts.py
index 1e0f317c..ff2b95da 100644
--- a/src/yaw/redshifts.py
+++ b/src/yaw/redshifts.py
@@ -137,7 +137,7 @@ def from_correlation_data(cls, *args, **kwargs):
         .. deprecated:: 2.3.2
             Renamed to :meth:`from_corrdata`.
         """
-        return cls.from_corrdata(*args, **kwargs)
+        return cls.from_corrdata(*args, **kwargs)  # pragma: no cover
 
     @classmethod
     def from_corrdata(
@@ -200,7 +200,7 @@ def from_corrdata(
             mitigate.append("unknown")
 
         if len(mitigate) > 0:
-            logger.debug(f"mitigating {' and '.join(mitigate)} sample bias")
+            logger.debug("mitigating %s sample bias", " and ".join(mitigate))
         N = cross_data.n_samples
         dzsq_data = cross_data.dz**2
         dzsq_samp = np.tile(dzsq_data, N).reshape((N, -1))
@@ -221,7 +221,7 @@ def from_correlation_functions(cls, *args, **kwargs):
         .. deprecated:: 2.3.2
             Renamed to :meth:`from_corrfuncs`.
         """
-        return cls.from_corrfuncs(*args, **kwargs)
+        return cls.from_corrfuncs(*args, **kwargs)  # pragma: no cover
 
     @classmethod
     def from_corrfuncs(
@@ -448,7 +448,7 @@ def shift(
         )
 
 
-@dataclass(frozen=True, repr=False, eq=False)
+@dataclass(frozen=True, repr=False)
 class HistData(RedshiftData):
     """Container for histogram data.
 
@@ -477,6 +477,12 @@ class HistData(RedshiftData):
     density: bool = field(default=False)
     """Whether the data is normalised, i.e. a density estimate."""
 
+    def __eq__(self, other: object) -> bool:
+        parent_eq = super().__eq__(other)
+        if isinstance(other, self.__class__):
+            return self.density == other.density
+        return parent_eq
+
     @property
     def _dat_desc(self) -> str:
         n = "normalised " if self.density else " "


Test patch:
diff --git a/tests/yaw/config/test_backend.py b/tests/yaw/config/test_backend.py
index 65390035..7c4eedf4 100644
--- a/tests/yaw/config/test_backend.py
+++ b/tests/yaw/config/test_backend.py
@@ -9,6 +9,16 @@ class TestBackendConfig:
     def test_default_threads(self):
         assert BackendConfig().thread_num == os.cpu_count()
 
+    def test_modify(self):
+        substitutes = dict(
+            thread_num=20,
+            crosspatch=False,
+            rbin_slop=-1.0,
+        )
+        for param, value in substitutes.items():
+            conf = BackendConfig().modify(**{param: value})
+            assert getattr(conf, param) == value
+
     def test_get_threads(self):
         conf = BackendConfig(thread_num=10)
         assert conf.get_threads() == 10
diff --git a/tests/yaw/config/test_binning.py b/tests/yaw/config/test_binning.py
index 6333194d..2e2d4594 100644
--- a/tests/yaw/config/test_binning.py
+++ b/tests/yaw/config/test_binning.py
@@ -1,12 +1,8 @@
 import numpy as np
 import numpy.testing as npt
-from pytest import fixture, raises
+from pytest import fixture, raises, warns
 
-from yaw.config.binning import (
-    AutoBinningConfig,
-    ManualBinningConfig,
-    make_binning_config,
-)
+from yaw.config.binning import BinningConfig
 from yaw.config.utils import ConfigError
 
 
@@ -28,64 +24,33 @@ def bin_props(bin_edges):
 
 @fixture
 def manual_binning(bin_edges):
-    return ManualBinningConfig(bin_edges)
+    return BinningConfig.create(zbins=bin_edges)
 
 
 @fixture
 def auto_binning(bin_props):
-    return AutoBinningConfig.generate(**bin_props)
+    return BinningConfig.create(**bin_props)
 
 
-class TestManualBinningConfig:
+class TestBinningConfig:
     def test_init(self, bin_edges, manual_binning):
         npt.assert_array_equal(bin_edges, manual_binning.zbins)
         assert manual_binning.method == "manual"
         with raises(ValueError):
-            ManualBinningConfig(bin_edges[::-1])
+            BinningConfig(zbins=bin_edges[::-1])
         with raises(ConfigError):
-            ManualBinningConfig([0.1])
+            BinningConfig(zbins=[0.1])
 
     def test_pass_repr(self, manual_binning):
         str(manual_binning)
 
-    def test_eq(self, bin_edges, manual_binning):
+    def test_eq(self, bin_edges, bin_props, auto_binning, manual_binning):
         new_bins = bin_edges.copy()
         new_bins[-1] = 2.0
-        assert ManualBinningConfig(bin_edges) == manual_binning
-        assert not ManualBinningConfig(new_bins) == manual_binning
+        assert BinningConfig.create(zbins=bin_edges) == manual_binning
+        assert not BinningConfig.create(zbins=new_bins) == manual_binning
 
-    def test_properties(self, bin_props, manual_binning):
-        assert manual_binning.zmin == bin_props["zmin"]
-        assert manual_binning.zmax == bin_props["zmax"]
-        assert manual_binning.zbin_num == bin_props["zbin_num"]
-
-    def test_to_dict(self, bin_edges, manual_binning):
-        the_dict = manual_binning.to_dict()
-        assert set(the_dict.keys()) == {"zbins", "method"}
-        npt.assert_array_equal(the_dict["zbins"], bin_edges)
-
-    def test_from_dict(self, bin_edges, manual_binning):
-        the_dict = dict(zbins=bin_edges)
-        assert ManualBinningConfig.from_dict(the_dict) == manual_binning
-        # provide arugments only consumed by AutoBinningConfig
-        the_dict.update(dict(zmin=bin_edges.min(), zmax=bin_edges.max()))
-        with raises(TypeError):
-            ManualBinningConfig.from_dict(the_dict)
-
-
-class TestAutoBinningConfig:
-    def test_generate(self, bin_edges, auto_binning):
-        method = "linear"
-        binning = AutoBinningConfig(bin_edges, method=method)
-        npt.assert_array_equal(bin_edges, binning.zbins)
-        assert binning.method == method
-        npt.assert_array_equal(bin_edges, auto_binning.zbins)
-
-    def test_pass_repr(self, auto_binning):
-        str(auto_binning)
-
-    def test_eq(self, bin_props, auto_binning):
-        assert AutoBinningConfig.generate(**bin_props) == auto_binning
+        assert BinningConfig.create(**bin_props) == auto_binning
         new_props = dict(
             zmin=bin_props["zmin"] + 0.001,
             zmax=bin_props["zmax"] + 0.001,
@@ -95,48 +60,82 @@ def test_eq(self, bin_props, auto_binning):
         for param, value in new_props.items():
             the_props = {k: v for k, v in bin_props.items()}
             the_props[param] = value
-            assert not AutoBinningConfig.generate(**the_props) == auto_binning
+            assert not BinningConfig.create(**the_props) == auto_binning
+
+        assert auto_binning != 1
+
+    def test_properties(self, bin_props, auto_binning, manual_binning):
+        assert manual_binning.zmin == bin_props["zmin"]
+        assert manual_binning.zmax == bin_props["zmax"]
+        assert manual_binning.zbin_num == bin_props["zbin_num"]
+        assert manual_binning.method == "manual"
+        assert manual_binning.is_manual
 
-    def test_properties(self, bin_props, auto_binning):
         assert auto_binning.zmin == bin_props["zmin"]
         assert auto_binning.zmax == bin_props["zmax"]
         assert auto_binning.zbin_num == bin_props["zbin_num"]
+        assert auto_binning.method == "linear"
+        assert not auto_binning.is_manual
 
-    def test_to_dict(self, auto_binning):
-        the_dict = auto_binning.to_dict()
-        assert set(the_dict.keys()) == {"zmin", "zmax", "zbin_num", "method"}
+    def test_create_auto(self, bin_props, auto_binning):
+        assert auto_binning == BinningConfig.create(**bin_props)
+        the_props = {k: v for k, v in bin_props.items()}
+        the_props["method"] = None
+        with raises(ValueError):
+            BinningConfig.create(**the_props)
 
-    def test_from_dict(self, bin_props, auto_binning):
-        the_dict = {key: val for key, val in bin_props.items()}
-        assert AutoBinningConfig.from_dict(the_dict) == auto_binning
-        # provide arugments only consumed by ManualBinningConfig
-        the_dict.update(dict(zbins=bin_edges))
-        with raises(TypeError):
-            AutoBinningConfig.from_dict(the_dict)
+    def test_create_manual(self, bin_edges, manual_binning):
+        assert manual_binning == BinningConfig.create(zbins=bin_edges)
 
+    def test_create_all_arg(self, bin_edges, bin_props, auto_binning):
+        with warns(UserWarning, match="ignored"):
+            assert BinningConfig.create(zbins=bin_edges, **bin_props) == auto_binning
 
-def test_make_binning_config_auto(bin_props, auto_binning):
-    assert auto_binning == make_binning_config(**bin_props)
-    the_props = {k: v for k, v in bin_props.items()}
-    the_props["method"] = None
-    with raises(ValueError):
-        make_binning_config(**the_props)
+    def test_create_no_z_arg(self, bin_edges, bin_props, auto_binning):
+        with raises(ConfigError):
+            BinningConfig.create(zmin=0.1)
+        with raises(ConfigError):
+            BinningConfig.create(zmax=0.1)
+        with raises(ConfigError):
+            BinningConfig.create()
 
+    def test_modify_auto(self, bin_props, auto_binning, manual_binning):
+        assert auto_binning == auto_binning.modify(**bin_props)
+        assert auto_binning == manual_binning.modify(**bin_props)
 
-def test_make_binning_config_manual(bin_edges, manual_binning):
-    assert manual_binning == make_binning_config(zbins=bin_edges)
+        for param, value in bin_props.items():
+            conf = auto_binning.modify(**{param: value})
+            assert getattr(conf, param) == value
 
+        for param in ("zmin", "zmax"):
+            with raises(ConfigError):
+                manual_binning.modify(**{param: bin_props[param]})
 
-def test_make_binning_config_all_arg(bin_edges, bin_props, manual_binning):
-    assert isinstance(
-        make_binning_config(zbins=bin_edges, **bin_props), AutoBinningConfig
-    )
+    def test_modify_manual(self, bin_edges, auto_binning, manual_binning):
+        assert manual_binning == auto_binning.modify(zbins=bin_edges)
+        assert manual_binning == manual_binning.modify(zbins=bin_edges)
+
+    def test_modify_all_arg(self, bin_edges, bin_props, auto_binning, manual_binning):
+        assert auto_binning.modify(zbins=bin_edges, **bin_props).is_manual
+        assert manual_binning.modify(zbins=bin_edges, **bin_props).is_manual
+
+    def test_modify_no_z_arg(self, manual_binning):
+        with raises(ConfigError):
+            manual_binning.modify(zmin=0.1)
+        with raises(ConfigError):
+            manual_binning.modify(zmax=0.1)
 
+    def test_to_dict(self, bin_edges, auto_binning, manual_binning):
+        the_dict = manual_binning.to_dict()
+        assert set(the_dict.keys()) == {"zbins", "method"}
+        npt.assert_array_equal(the_dict["zbins"], bin_edges)
 
-def test_make_binning_config_no_z_arg(bin_edges, bin_props, auto_binning):
-    with raises(ConfigError):
-        make_binning_config(zmin=0.1)
-    with raises(ConfigError):
-        make_binning_config(zmax=0.1)
-    with raises(ConfigError):
-        make_binning_config()
+        the_dict = auto_binning.to_dict()
+        assert set(the_dict.keys()) == {"zmin", "zmax", "zbin_num", "method"}
+
+    def test_from_dict(self, bin_edges, bin_props, auto_binning, manual_binning):
+        the_dict = dict(zbins=bin_edges)
+        assert BinningConfig.from_dict(the_dict) == manual_binning
+
+        the_dict = {key: val for key, val in bin_props.items()}
+        assert BinningConfig.from_dict(the_dict) == auto_binning
diff --git a/tests/yaw/config/test_configuration.py b/tests/yaw/config/test_configuration.py
index 3f9d4d17..3f5c30c1 100644
--- a/tests/yaw/config/test_configuration.py
+++ b/tests/yaw/config/test_configuration.py
@@ -3,12 +3,7 @@
 from astropy.cosmology import WMAP9, Planck15
 from pytest import fixture, mark, raises
 
-from yaw.config import (
-    AutoBinningConfig,
-    BackendConfig,
-    Configuration,
-    ManualBinningConfig,
-)
+from yaw.config import BackendConfig, BinningConfig, Configuration
 from yaw.config.utils import ConfigError
 from yaw.core.cosmology import get_default_cosmology
 
@@ -185,27 +180,29 @@ def test_update_binning_auto_update_other(self, default_config):
             conf = auto_config.modify(**{param: value})
             auto_kwargs = auto_config.to_dict()["binning"]
             auto_kwargs[param] = value
-            assert conf.binning == AutoBinningConfig.generate(**auto_kwargs)
+            assert conf.binning == BinningConfig.create(**auto_kwargs)
 
     def test_update_binning_manual_update_zbins(self, manual_config):
         zbins = np.linspace(1, 10)
-        assert manual_config.modify(zbins=zbins).binning == ManualBinningConfig(zbins)
+        assert manual_config.modify(zbins=zbins).binning == BinningConfig.create(
+            zbins=zbins
+        )
 
     def test_update_binning_manual_all_other(self, manual_config, default_dict):
         bin_dict = default_dict["binning"]
         kwargs = dict(zmin=bin_dict["zmin"], zmax=bin_dict["zmax"])
         # only z limits
-        expect = AutoBinningConfig.generate(**bin_dict)
+        expect = BinningConfig.create(**bin_dict)
         got = manual_config.modify(**kwargs)
         assert expect == got.binning
         # add method
         added_dict = {k: v for k, v in bin_dict.items() if k != "method"}
-        expect = AutoBinningConfig.generate(**added_dict, method="comoving")
+        expect = BinningConfig.create(**added_dict, method="comoving")
         got = manual_config.modify(**kwargs, method="comoving")
         assert expect == got.binning
         # add zbin_num
         added_dict = {k: v for k, v in bin_dict.items() if k != "zbin_num"}
-        expect = AutoBinningConfig.generate(**added_dict, zbin_num=5)
+        expect = BinningConfig.create(**added_dict, zbin_num=5)
         got = manual_config.modify(**kwargs, zbin_num=5)
         assert expect == got.binning
 
diff --git a/tests/yaw/config/test_resampling.py b/tests/yaw/config/test_resampling.py
index 6e77df4b..68904454 100644
--- a/tests/yaw/config/test_resampling.py
+++ b/tests/yaw/config/test_resampling.py
@@ -19,6 +19,18 @@ def test_n_patches(self, method):
         conf.get_samples(n_patches)
         assert conf.n_patches == n_patches
 
+    def test_modify(self):
+        substitutes = dict(
+            method="bootstrap",
+            crosspatch=False,
+            n_boot=10,
+            global_norm=True,
+            seed=99,
+        )
+        for param, value in substitutes.items():
+            conf = ResamplingConfig().modify(**{param: value})
+            assert getattr(conf, param) == value
+
     def test__generate_jackknife(self):
         samples = ResamplingConfig()._generate_jackknife(n_patches=3)
         npt.assert_array_equal(samples.sum(axis=1), [3, 2, 1])
diff --git a/tests/yaw/config/test_scales.py b/tests/yaw/config/test_scales.py
index bb948828..3293d997 100644
--- a/tests/yaw/config/test_scales.py
+++ b/tests/yaw/config/test_scales.py
@@ -31,9 +31,12 @@ def test_init_mixed(self):
         with raises(ConfigError):
             ScalesConfig(100, [100])
 
-    def test_dict_keys(self, default_scales):
-        scales = set(default_scales.dict_keys())
-        assert scales == set(str(scale) for scale in default_scales)
+    def test_iter(self, default_scales):
+        scales = set(
+            Scale(rmin, rmax)
+            for rmin, rmax in zip(default_scales.rmin, default_scales.rmax)
+        )
+        assert scales == set(default_scales)
 
     def test_getitem(self, default_scales):
         assert default_scales[1] == Scale(
@@ -45,3 +48,15 @@ def test_eq(self, default_scales):
         assert scales != ScalesConfig(scales.rmin[:1], scales.rmax[:1])
         assert scales != ScalesConfig(scales.rmin, scales.rmax, rweight=1.0)
         assert scales != ScalesConfig(scales.rmin, scales.rmax, rbin_num=11)
+        assert scales != 1
+
+    def test_modify(self, default_scales):
+        substitutes = dict(
+            rmin=[10, 20],
+            rmax=[200, 400],
+            rweight=-1.0,
+            rbin_num=100,
+        )
+        for param, value in substitutes.items():
+            conf = default_scales.modify(**{param: value})
+            assert getattr(conf, param) == value
diff --git a/tests/yaw/core/test_containers.py b/tests/yaw/core/test_containers.py
new file mode 100644
index 00000000..5d1715d0
--- /dev/null
+++ b/tests/yaw/core/test_containers.py
@@ -0,0 +1,76 @@
+import numpy as np
+import numpy.testing as npt
+import pandas as pd
+import pandas.testing as pdt
+from pytest import fixture, raises
+
+from yaw.core import containers
+
+
+@fixture
+def binning():
+    return pd.IntervalIndex.from_breaks([0.0, 0.5, 1.0, 1.5])
+
+
+@fixture
+def test_data_samples(binning):
+    n_bins, n_samples = len(binning), 7
+    data = np.full(n_bins, 4)
+    samples = np.repeat(np.arange(1, 8), n_bins).reshape((n_samples, n_bins))
+    return binning, data, samples
+
+
+class TestSampledData:
+    def test_init(self, test_data_samples):
+        binning, data, samples = test_data_samples
+        sd = containers.SampledData(binning, data, samples, method="jackknife")
+        # binning and data mismatch
+        with raises(ValueError, match="unexpected *"):
+            containers.SampledData(binning, data[1:], samples, method="jackknife")
+        # wrong sample shape
+        with raises(ValueError, match="number of bins*"):
+            samples1 = samples[:, :1]
+            containers.SampledData(binning, data, samples1, method="jackknife")
+        # check some of the representation
+        assert f"n_bins={samples.shape[1]}" in repr(sd)
+        assert f"n_samples={samples.shape[0]}" in repr(sd)  # tests .n_samples
+        # wrong sampling method
+        with raises(ValueError):
+            containers.SampledData(binning, data, samples, method="garbage")
+
+    def test_getters(self, test_data_samples):
+        binning, data, samples = test_data_samples
+        sd = containers.SampledData(binning, data, samples, method="jackknife")
+        # test shape and index
+        pdt.assert_index_equal(binning, sd.get_data().index)
+        pdt.assert_index_equal(binning, sd.get_samples().index)
+        assert sd.get_data().shape == data.shape
+        assert sd.get_samples().shape == samples.T.shape
+
+    def test_compatible(self, test_data_samples):
+        binning, data, samples = test_data_samples
+        sd = containers.SampledData(binning, data, samples, method="jackknife")
+        # trivial case
+        assert sd.is_compatible(sd)
+        # wrong type
+        with raises(TypeError):
+            sd.is_compatible(1)
+        # different binning
+        binning2 = pd.IntervalIndex.from_breaks(np.linspace(1.0, 2.0, len(binning) + 1))
+        sd2 = containers.SampledData(binning2, data, samples, method="jackknife")
+        assert not sd.is_compatible(sd2)
+        # different number of samples
+        sd2 = containers.SampledData(binning, data, samples[:-1], method="jackknife")
+        assert not sd.is_compatible(sd2)
+        # different method
+        sd2 = containers.SampledData(binning, data, samples, method="bootstrap")
+        assert not sd.is_compatible(sd2)
+
+    def test_binning(self, test_data_samples):
+        # methods defined in BinnedQuantity
+        binning, data, samples = test_data_samples
+        sd = containers.SampledData(binning, data, samples, method="jackknife")
+        assert sd.n_bins == len(data)
+        npt.assert_equal(sd.mids, np.array([0.25, 0.75, 1.25]))
+        npt.assert_equal(sd.edges, np.array([0.0, 0.5, 1.0, 1.5]))
+        npt.assert_equal(sd.dz, np.array([0.5, 0.5, 0.5]))
diff --git a/tests/yaw/core/test_coordinates.py b/tests/yaw/core/test_coordinates.py
index 90f8dcbf..92ab6fca 100644
--- a/tests/yaw/core/test_coordinates.py
+++ b/tests/yaw/core/test_coordinates.py
@@ -310,6 +310,10 @@ def test_iter(self):
         for a, b in zip(dists, dists.values):
             assert a == coordinates.Dist3D(b)
 
+    def test_eq(self):
+        assert coordinates.Dist3D(1.0) == coordinates.Dist3D(1.0)
+        assert coordinates.Dist3D(1.0) != 1.0
+
     def test_ordering(self):
         assert coordinates.Dist3D(1.0) == coordinates.Dist3D(1.0)
         assert coordinates.Dist3D(1.0) <= coordinates.Dist3D(1.0)
@@ -319,6 +323,8 @@ def test_ordering(self):
         assert coordinates.Dist3D(2.0) > coordinates.Dist3D(1.0)
         # test array
         assert np.all(coordinates.Dist3D([0.0, 1.0]) == coordinates.Dist3D([0.0, 1.0]))
+        with raises(TypeError):
+            coordinates.Dist3D(1.0) < 1.0
 
     def test_add_sub(self):
         dist_90deg = coordinates.Dist3D(np.sqrt(2.0))  # don't add linearly
@@ -335,6 +341,11 @@ def test_add_sub(self):
         )
         # wrap around unit sphere
         assert coordinates.Dist3D(2.0) + dist_90deg == dist_90deg
+        # test incompatible types
+        with raises(TypeError):
+            coordinates.Dist3D(1.0) + 1.0
+        with raises(TypeError):
+            coordinates.Dist3D(1.0) - 1.0
 
     def test_conversion(self):
         assert coordinates.Dist3D(0.0).to_sky() == coordinates.DistSky(0.0)
@@ -374,6 +385,13 @@ def test_values(self):
         dists.__repr__()
         f"{coordinates.DistSky(2.0)}"
 
+    def test_add_sub(self):
+        # test incompatible types
+        with raises(TypeError):
+            coordinates.DistSky(1.0) + 1.0
+        with raises(TypeError):
+            coordinates.DistSky(1.0) - 1.0
+
     def test_conversion(self):
         # see also TestDist3D.test_to_sky
         assert coordinates.DistSky(0.0).to_3d() == coordinates.Dist3D(0.0)
diff --git a/tests/yaw/core/test_cosmology.py b/tests/yaw/core/test_cosmology.py
index 796257d0..11b742b2 100644
--- a/tests/yaw/core/test_cosmology.py
+++ b/tests/yaw/core/test_cosmology.py
@@ -25,6 +25,9 @@ def test_init(self):
         with raises(ValueError):
             cosmology.Scale(10, 2)
 
+    def test_str(self):
+        str(cosmology.Scale(1.0, 100.0))
+
     def test_mid(self):
         assert cosmology.Scale(10.0, 30.0).mid == 20.0
 
diff --git a/tests/yaw/correlation/test_corrfunc.py b/tests/yaw/correlation/test_corrfunc.py
index 98d57799..d2466d38 100644
--- a/tests/yaw/correlation/test_corrfunc.py
+++ b/tests/yaw/correlation/test_corrfunc.py
@@ -1,9 +1,10 @@
 import numpy as np
 import numpy.testing as npt
 import pandas as pd
-from pytest import fixture
+from pytest import fixture, mark, raises
 
-from yaw.correlation import corrfuncs
+from yaw import examples
+from yaw.correlation import corrfuncs, estimators
 
 
 @fixture
@@ -45,7 +46,15 @@ def stats_bootstrap(binning):
 @fixture
 def correlation_data(test_data_samples):
     binning, data, samples = test_data_samples
-    return corrfuncs.CorrData(binning, data, samples, method="jackknife")
+    return corrfuncs.CorrData(binning, data, samples, method="jackknife", info="info")
+
+
+def text_grep(path, pattern, sub):
+    with open(path) as f:
+        lines = [line.replace(pattern, sub) for line in f.readlines()]
+    with open(path, "w") as f:
+        for line in lines:
+            f.write(line)
 
 
 class TestCorrData:
@@ -70,3 +79,261 @@ def test_io(self, tmp_path_factory, correlation_data):
         restored = correlation_data.from_files(path)
         npt.assert_allclose(correlation_data.data, restored.data)
         npt.assert_allclose(correlation_data.samples, restored.samples)
+        assert correlation_data.info == restored.info
+        # break the data file
+        text_grep(f"{path}.smp", correlation_data.method[:4], "udef")
+        with raises(ValueError):
+            correlation_data.from_files(path)
+        text_grep(f"{path}.smp", "udef", correlation_data.method[:4])
+        text_grep(f"{path}.smp", "z_low", "wrong_key")
+        with raises(ValueError):
+            correlation_data.from_files(path)
+
+    @mark.parametrize("flags", [False, True])
+    def test_plot(self, correlation_data, flags):
+        correlation_data.plot(error_bars=flags, zero_line=flags, scale_by_dz=flags)
+
+    @mark.parametrize("redshift", [False, True])
+    def test_plot_corr(self, correlation_data, redshift):
+        correlation_data.plot_corr(redshift=redshift)
+
+
+def test_check_mergable():
+    corrfuncs.check_mergable([examples.w_sp, examples.w_sp])
+    # different pair counts
+    with raises(ValueError, match="incompatible"):
+        corrfuncs.check_mergable([examples.w_sp, examples.w_ss])
+
+
+class TestCorrFunc:
+    def test_init(self):
+        pairs = dict(
+            dd=examples.w_ss.dd,
+            dr=examples.w_ss.dr,
+            rd=examples.w_ss.dr,  # duplicated
+            rr=examples.w_ss.rr,
+        )
+        assert examples.w_ss == corrfuncs.CorrFunc(
+            dd=pairs["dd"], dr=pairs["dr"], rr=pairs["rr"]
+        )
+        with raises(ValueError, match="required"):
+            corrfuncs.CorrFunc(dd=pairs["dd"])
+        for name, pair in pairs.items():
+            if name == "dd":
+                continue
+            corrfuncs.CorrFunc(dd=pairs["dd"], **{name: pair})
+        with raises(ValueError):
+            corrfuncs.CorrFunc(dd=pairs["dd"], dr=pairs["dr"].bins[1:])
+        with raises(ValueError):
+            corrfuncs.CorrFunc(dd=pairs["dd"], dr=pairs["dr"].patches[1:])
+
+    def test_repr(self):
+        repr(examples.w_ss)
+
+    def test_eq(self):
+        assert examples.w_sp == examples.w_sp
+        assert examples.w_sp != examples.w_ss
+        assert examples.w_sp != 1
+
+    def test_add(self):
+        summed = examples.w_ss + examples.w_ss
+        assert summed.dd == examples.w_ss.dd + examples.w_ss.dd
+        assert summed.dr == examples.w_ss.dr + examples.w_ss.dr
+        assert summed.rr == examples.w_ss.rr + examples.w_ss.rr
+        with raises(ValueError, match="not set"):
+            new = corrfuncs.CorrFunc(dd=summed.dd, dr=summed.dr)
+            summed + new
+        with raises(TypeError):
+            summed + 1
+
+    def test_radd(self):
+        assert (0 + examples.w_sp) == examples.w_sp
+        assert sum([examples.w_sp]) == examples.w_sp
+        with raises(TypeError):
+            1 + examples.w_sp
+
+    def test_mul(self):
+        mul = examples.w_sp * 2
+        assert mul == examples.w_sp + examples.w_sp
+        assert mul.dd.total == examples.w_sp.dd.total
+        assert mul.dd.count == examples.w_sp.dd.count * 2
+        with raises(TypeError):
+            examples.w_sp * True
+
+    def test_auto(self):
+        assert examples.w_sp.auto == examples.w_sp.dd.count.auto
+
+    def test_n_patches(self):
+        assert examples.w_sp.n_patches == examples.w_sp.dd.count.n_patches
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_bins(self, items):
+        subset = examples.w_sp.bins[items]
+        assert subset == corrfuncs.CorrFunc(
+            dd=examples.w_sp.dd.bins[items], dr=examples.w_sp.dr.bins[items]
+        )
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_patches(self, items):
+        subset = examples.w_sp.patches[items]
+        assert subset == corrfuncs.CorrFunc(
+            dd=examples.w_sp.dd.patches[items], dr=examples.w_sp.dr.patches[items]
+        )
+
+    def test_get_binning(self):
+        assert (examples.w_sp.get_binning() == examples.w_sp.dd.get_binning()).all()
+
+    def test_is_compatible(self):
+        examples.w_sp.is_compatible(examples.w_sp)
+        with raises(ValueError):
+            assert not examples.w_sp.is_compatible(examples.w_sp.bins[1:])
+            assert examples.w_sp.is_compatible(examples.w_sp.bins[1:], require=True)
+        with raises(ValueError):
+            assert not examples.w_sp.is_compatible(examples.w_sp.patches[1:])
+            assert examples.w_sp.is_compatible(examples.w_sp.patches[1:], require=True)
+
+    def test_estimators(self):
+        pairs = dict(
+            dd=examples.w_ss.dd,
+            dr=examples.w_ss.dr,
+            rd=examples.w_ss.dr,  # duplicated
+            rr=examples.w_ss.rr,
+        )
+        cf = corrfuncs.CorrFunc(**pairs).estimators
+        assert set(cf.keys()) == {"HM", "LS", "DP", "PH"}
+        cf = corrfuncs.CorrFunc(
+            dd=pairs["dd"], dr=pairs["dr"], rr=pairs["rr"]
+        ).estimators
+        assert set(cf.keys()) == {"HM", "LS", "DP", "PH"}
+        cf = corrfuncs.CorrFunc(
+            dd=pairs["dd"], rd=pairs["rd"], rr=pairs["rr"]
+        ).estimators
+        assert set(cf.keys()) == {"DP", "PH"}
+
+        cf = corrfuncs.CorrFunc(
+            dd=pairs["dd"], dr=pairs["dr"], rd=pairs["rd"]
+        ).estimators
+        assert set(cf.keys()) == {"DP"}
+        cf = corrfuncs.CorrFunc(dd=pairs["dd"], rd=pairs["rd"]).estimators
+        assert set(cf.keys()) == {"DP"}
+        cf = corrfuncs.CorrFunc(dd=pairs["dd"], dr=pairs["dr"]).estimators
+        assert set(cf.keys()) == {"DP"}
+
+        cf = corrfuncs.CorrFunc(dd=pairs["dd"], rr=pairs["rr"]).estimators
+        assert set(cf.keys()) == {"PH"}
+
+    def test__check_and_select_estimator(self):
+        pairs = dict(
+            dd=examples.w_ss.dd,
+            dr=examples.w_ss.dr,
+            rd=examples.w_ss.dr,  # duplicated
+            rr=examples.w_ss.rr,
+        )
+
+        cf = corrfuncs.CorrFunc(**pairs)
+        with raises(ValueError, match="invalid"):
+            cf._check_and_select_estimator("undefined")
+        assert cf._check_and_select_estimator() == estimators.LandySzalay
+        for est in (
+            estimators.LandySzalay,
+            estimators.Hamilton,
+            estimators.DavisPeebles,
+            estimators.PeeblesHauser,
+        ):
+            assert cf._check_and_select_estimator(est.short) == est
+
+        for pair in ("dr", "rd"):
+            cf = corrfuncs.CorrFunc(dd=pairs["dd"], **{pair: pairs[pair]})
+            assert cf._check_and_select_estimator() == estimators.DavisPeebles
+            for est in (estimators.DavisPeebles,):
+                assert cf._check_and_select_estimator(est.short) == est
+            for est in (
+                estimators.LandySzalay,
+                estimators.Hamilton,
+                estimators.PeeblesHauser,
+            ):
+                with raises(estimators.EstimatorError, match="requires"):
+                    cf._check_and_select_estimator(est.short)
+
+        cf = corrfuncs.CorrFunc(dd=pairs["dd"], rr=pairs["rr"])
+        assert cf._check_and_select_estimator() == estimators.PeeblesHauser
+        for est in (estimators.PeeblesHauser,):
+            assert cf._check_and_select_estimator(est.short) == est
+        for est in (
+            estimators.LandySzalay,
+            estimators.Hamilton,
+            estimators.DavisPeebles,
+        ):
+            with raises(estimators.EstimatorError, match="requires"):
+                cf._check_and_select_estimator(est.short)
+
+    def test__getattr_from_cts(self):
+        pairs = dict(
+            dd=examples.w_ss.dd,
+            dr=examples.w_ss.dr,
+            rd=examples.w_ss.dr,  # duplicated
+            rr=examples.w_ss.rr,
+        )
+        with_rd = corrfuncs.CorrFunc(**pairs)
+        rd = pairs.pop("rd")
+        without_rd = corrfuncs.CorrFunc(**pairs)
+        # CtsDD
+        assert with_rd.dd == with_rd._getattr_from_cts(estimators.CtsDD())
+        assert with_rd.dd == without_rd._getattr_from_cts(estimators.CtsDD())
+        # CtsDR
+        assert with_rd.dr == with_rd._getattr_from_cts(estimators.CtsDR())
+        assert with_rd.dr == without_rd._getattr_from_cts(estimators.CtsDR())
+        # CtsRD
+        assert with_rd.rd == with_rd._getattr_from_cts(estimators.CtsRD())
+        assert without_rd._getattr_from_cts(estimators.CtsRD()) is None
+        # CtsRR
+        assert with_rd.rr == with_rd._getattr_from_cts(estimators.CtsRR())
+        assert with_rd.rr == without_rd._getattr_from_cts(estimators.CtsRR())
+        # CtsMix
+        assert with_rd.dr == with_rd._getattr_from_cts(estimators.CtsMix())
+        assert with_rd.dr == without_rd._getattr_from_cts(estimators.CtsMix())
+        without_dr = corrfuncs.CorrFunc(**pairs, rd=rd)
+        assert without_dr.rd == without_dr._getattr_from_cts(estimators.CtsMix())
+
+    @mark.skip
+    def test_sample(self):
+        assert 0
+
+    @mark.parametrize("cf", [examples.w_sp, examples.w_ss, examples.w_pp])
+    def test_io(self, tmp_path_factory, cf):
+        path = tmp_path_factory.mktemp("corr_func_file") / "testfile.hdf"
+        cf.to_file(path)
+        assert corrfuncs.CorrFunc.from_file(path) == cf
+
+    @mark.skip
+    def test_concatenate_patches(self):
+        assert 0
+
+    @mark.skip
+    def test_concatenate_bins(self):
+        assert 0
+
+
+@mark.skip
+def test__create_dummy_counts(self):
+    assert 0
+
+
+@mark.skip
+def test_add_corrfuncs(self):
+    assert 0
+
+
+@mark.skip
+def test__check_patch_centers(self):
+    assert 0
+
+
+@mark.skip
+def test_autocorrelate(self):
+    assert 0
+
+
+@mark.skip
+def test_crosscorrelate(self):
+    assert 0
diff --git a/tests/yaw/correlation/test_estimators.py b/tests/yaw/correlation/test_estimators.py
new file mode 100644
index 00000000..e93a0fc6
--- /dev/null
+++ b/tests/yaw/correlation/test_estimators.py
@@ -0,0 +1,67 @@
+import numpy as np
+from pytest import mark, raises, warns
+
+from yaw.correlation import estimators
+
+
+class TestCts:
+    def test_set_squashing(self):
+        assert len({estimators.CtsDR(), estimators.CtsRD(), estimators.CtsMix()}) == 2
+
+    def test_equivalents(self):
+        assert estimators.CtsDR() == estimators.CtsMix()
+        assert estimators.CtsRD() == estimators.CtsMix()
+        assert estimators.CtsDD() != 1.0
+
+    def test_repr_checks(self):
+        repr(estimators.CtsMix())
+        str(estimators.CtsMix())
+
+
+@mark.parametrize(
+    "code,counts",
+    [
+        ("dd", estimators.CtsDD),
+        ("dr", estimators.CtsDR),
+        ("rd", estimators.CtsRD),
+        ("rr", estimators.CtsRR),
+    ],
+)
+def test_cts_from_code(code, counts):
+    assert counts() == estimators.cts_from_code(code)
+    with raises(ValueError):
+        estimators.cts_from_code("undefined symbol")
+
+
+class TestPeeblesHauser:
+    def test_eval(self):
+        assert estimators.PeeblesHauser.eval(dd=1, rr=1) == 0
+        with warns(UserWarning):
+            estimators.PeeblesHauser.eval(dd=np.float_(1), rr=np.float_(0))
+
+
+class TestDavisPeebles:
+    def test_eval(self):
+        assert estimators.DavisPeebles.eval(dd=1, dr_rd=1) == 0
+        with warns(UserWarning):
+            estimators.DavisPeebles.eval(dd=np.float_(1), dr_rd=np.float_(0))
+
+
+class TestHamilton:
+    def test_eval(self):
+        assert estimators.Hamilton.eval(dd=1, dr=1, rr=1) == 0
+        assert estimators.Hamilton.eval(dd=1, dr=1, rd=1, rr=1) == 0
+        with warns(UserWarning):
+            estimators.PeeblesHauser.eval(
+                dd=np.float_(1), dr=np.float_(1), rr=np.float_(0)
+            )
+
+
+class TestLandySzalay:
+    def test_eval(self):
+        assert estimators.LandySzalay.eval(dd=1, dr=1, rr=1) == 0
+        assert estimators.LandySzalay.eval(dd=1, dr=1, rd=1, rr=1) == 0
+        with warns(UserWarning):
+            estimators.PeeblesHauser.eval(
+                dd=np.float_(1), dr=np.float_(1), rr=np.float_(0)
+            )
diff --git a/tests/yaw/correlation/test_paircounts.py b/tests/yaw/correlation/test_paircounts.py
index 6b08a09a..1da586c1 100644
--- a/tests/yaw/correlation/test_paircounts.py
+++ b/tests/yaw/correlation/test_paircounts.py
@@ -8,9 +8,10 @@
 import numpy.testing as npt
 import pandas as pd
 import pandas.testing as pdt
-from pytest import fixture, raises
+from pytest import fixture, mark, raises
 
 from yaw.config import ResamplingConfig
+from yaw.core.math import apply_slice_ndim
 from yaw.correlation import paircounts
 
 
@@ -136,67 +137,63 @@ def binning():
 
 
 @fixture
-def test_data_samples(binning):
-    n_bins, n_samples = len(binning), 7
-    data = np.full(n_bins, 4)
-    samples = np.repeat(np.arange(1, 8), n_bins).reshape((n_samples, n_bins))
-    return binning, data, samples
-
-
-class TestSampledData:
-    def test_init(self, test_data_samples):
-        binning, data, samples = test_data_samples
-        sd = paircounts.SampledData(binning, data, samples, method="jackknife")
-        # binning and data mismatch
-        with raises(ValueError, match="unexpected *"):
-            paircounts.SampledData(binning, data[1:], samples, method="jackknife")
-        # wrong sample shape
-        with raises(ValueError, match="number of bins*"):
-            samples1 = samples[:, :1]
-            paircounts.SampledData(binning, data, samples1, method="jackknife")
-        # check some of the representation
-        assert f"n_bins={samples.shape[1]}" in repr(sd)
-        assert f"n_samples={samples.shape[0]}" in repr(sd)  # tests .n_samples
-        # wrong sampling method
-        with raises(ValueError):
-            paircounts.SampledData(binning, data, samples, method="garbage")
-
-    def test_getters(self, test_data_samples):
-        binning, data, samples = test_data_samples
-        sd = paircounts.SampledData(binning, data, samples, method="jackknife")
-        # test shape and index
-        pdt.assert_index_equal(binning, sd.get_data().index)
-        pdt.assert_index_equal(binning, sd.get_samples().index)
-        assert sd.get_data().shape == data.shape
-        assert sd.get_samples().shape == samples.T.shape
-
-    def test_compatible(self, test_data_samples):
-        binning, data, samples = test_data_samples
-        sd = paircounts.SampledData(binning, data, samples, method="jackknife")
-        # trivial case
-        assert sd.is_compatible(sd)
-        # wrong type
-        with raises(TypeError):
-            sd.is_compatible(1)
-        # different binning
-        binning2 = pd.IntervalIndex.from_breaks(np.linspace(1.0, 2.0, len(binning) + 1))
-        sd2 = paircounts.SampledData(binning2, data, samples, method="jackknife")
-        assert not sd.is_compatible(sd2)
-        # different number of samples
-        sd2 = paircounts.SampledData(binning, data, samples[:-1], method="jackknife")
-        assert not sd.is_compatible(sd2)
-        # different method
-        sd2 = paircounts.SampledData(binning, data, samples, method="bootstrap")
-        assert not sd.is_compatible(sd2)
-
-    def test_binning(self, test_data_samples):
-        # methods defined in BinnedQuantity
-        binning, data, samples = test_data_samples
-        sd = paircounts.SampledData(binning, data, samples, method="jackknife")
-        assert sd.n_bins == len(data)
-        npt.assert_equal(sd.mids, np.array([0.25, 0.75, 1.25]))
-        npt.assert_equal(sd.edges, np.array([0.0, 0.5, 1.0, 1.5]))
-        npt.assert_equal(sd.dz, np.array([0.5, 0.5, 0.5]))
+def binning_next():
+    return pd.IntervalIndex.from_breaks([1.5, 2.0, 2.5, 3.0])
+
+
+def modify_totals(totals, **kwargs):
+    allargs = dict(
+        binning=totals.get_binning(),
+        totals1=totals.totals1,
+        totals2=totals.totals2,
+        auto=totals.auto,
+    )
+    allargs.update(kwargs)
+    return paircounts.PatchedTotal(**allargs)
+
+
+def modify_counts(counts, **kwargs):
+    allargs = dict(
+        binning=counts.get_binning(),
+        counts=counts.counts,
+        auto=counts.auto,
+    )
+    allargs.update(kwargs)
+    return paircounts.PatchedCount(**allargs)
+
+
+def test_check_mergable_bins(patched_totals_full):
+    paircounts.check_mergable([patched_totals_full, patched_totals_full], patches=False)
+    new = modify_totals(patched_totals_full, auto=True)
+    with raises(ValueError, match="cross-"):
+        paircounts.check_mergable([patched_totals_full, new], patches=False)
+    with raises(ValueError, match="patch numbers"):
+        paircounts.check_mergable(
+            [patched_totals_full, patched_totals_full.patches[1:]], patches=False
+        )
+    paircounts.check_mergable(
+        [patched_totals_full, patched_totals_full.bins[1:]], patches=False
+    )
+
+
+def test_check_mergable_patches(patched_totals_full):
+    paircounts.check_mergable([patched_totals_full, patched_totals_full], patches=True)
+    new = modify_totals(patched_totals_full, auto=True)
+    with raises(ValueError, match="cross-"):
+        paircounts.check_mergable([patched_totals_full, new], patches=False)
+    paircounts.check_mergable(
+        [patched_totals_full, patched_totals_full.bins[1:]], patches=False
+    )
+
+
+def test_patch_idx_offset(patched_totals_full):
+    npt.assert_array_almost_equal(
+        paircounts.patch_idx_offset([patched_totals_full]), np.zeros(1)
+    )
+    result = np.array([0, patched_totals_full.n_patches])
+    npt.assert_array_almost_equal(
+        paircounts.patch_idx_offset([patched_totals_full, patched_totals_full]), result
+    )
 
 
 def test_binning_hdf(tmp_hdf5):
@@ -263,6 +260,29 @@ def test_init(self, binning, patch_totals, patched_totals_full):
         # just call once
         repr(pt)
 
+    def test_eq(self, patched_totals_full):
+        assert patched_totals_full == patched_totals_full
+        patched_new = paircounts.PatchedTotal(
+            patched_totals_full.get_binning(),
+            totals1=patched_totals_full.totals2,
+            totals2=patched_totals_full.totals2,
+            auto=patched_totals_full.auto,
+        )
+        assert patched_totals_full != patched_new
+        patched_new = paircounts.PatchedTotal(
+            patched_totals_full.get_binning(),
+            totals1=patched_totals_full.totals1,
+            totals2=patched_totals_full.totals2,
+            auto=not patched_totals_full.auto,
+        )
+        assert patched_totals_full != patched_new
+        assert patched_totals_full != 1
+
+    def test_size(self, patched_totals_full):
+        np = patched_totals_full.n_patches
+        nb = patched_totals_full.n_bins
+        assert patched_totals_full.size == (np * np * nb)
+
     def test_array(self, patched_totals_full, patch_matrix_full):
         pt = patched_totals_full
         n_bins = pt.n_bins
@@ -296,10 +316,65 @@ def test_jackknife(
             assert result.data[0] == data
             npt.assert_equal(result.samples[:, 0], jack)
 
+        subset = patched_totals_full.patches[0]
+        assert subset.sample_sum(config).samples.ndim == 2
+
     def test_bootstrap(self, patched_totals_full):
         with raises(NotImplementedError):
             patched_totals_full.sample_sum(ResamplingConfig(method="bootstrap"))
 
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_bins(self, patched_totals_full, patch_totals, items):
+        t1 = add_zbin_dimension(patch_totals[0], patched_totals_full.n_bins)
+        t2 = add_zbin_dimension(patch_totals[1], patched_totals_full.n_bins)
+        subset = patched_totals_full.bins[items]
+        npt.assert_array_equal(subset.totals1, apply_slice_ndim(t1, items, axis=1))
+        npt.assert_array_equal(subset.totals2, apply_slice_ndim(t2, items, axis=1))
+        assert (subset.get_binning() == patched_totals_full.get_binning()[items]).all()
+
+    def test_bins_iter(self, patched_totals_full, patch_totals):
+        t1 = add_zbin_dimension(patch_totals[0], patched_totals_full.n_bins)
+        t2 = add_zbin_dimension(patch_totals[1], patched_totals_full.n_bins)
+        for i, _bin in enumerate(patched_totals_full.bins):
+            npt.assert_array_equal(_bin.totals1, apply_slice_ndim(t1, i, axis=1))
+            npt.assert_array_equal(_bin.totals2, apply_slice_ndim(t2, i, axis=1))
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_patches(self, patched_totals_full, patch_totals, items):
+        t1 = add_zbin_dimension(patch_totals[0], patched_totals_full.n_bins)
+        t2 = add_zbin_dimension(patch_totals[1], patched_totals_full.n_bins)
+        subset = patched_totals_full.patches[items]
+        npt.assert_array_equal(subset.totals1, apply_slice_ndim(t1, items, axis=0))
+        npt.assert_array_equal(subset.totals2, apply_slice_ndim(t2, items, axis=0))
+        assert (subset.get_binning() == patched_totals_full.get_binning()).all()
+
+    def test_patches_iter(self, patched_totals_full, patch_totals):
+        t1 = add_zbin_dimension(patch_totals[0], patched_totals_full.n_bins)
+        t2 = add_zbin_dimension(patch_totals[1], patched_totals_full.n_bins)
+        for i, _bin in enumerate(patched_totals_full.patches):
+            npt.assert_array_equal(_bin.totals1, apply_slice_ndim(t1, i, axis=0))
+            npt.assert_array_equal(_bin.totals2, apply_slice_ndim(t2, i, axis=0))
+
+    def test_concatenate_patches(self, patched_totals_full, binning_next):
+        merged = patched_totals_full.concatenate_patches(patched_totals_full)
+        diag = np.diag(patched_totals_full.as_array()[:, :, 0])
+        npt.assert_array_equal(
+            np.diag(merged.as_array()[:, :, 0]), np.concatenate([diag, diag])
+        )
+        mod = modify_totals(patched_totals_full, binning=binning_next)
+        with raises(ValueError):
+            patched_totals_full.concatenate_patches(mod)
+
+    def test_concatenate_bins(self, patched_totals_full, binning_next):
+        appended = modify_totals(patched_totals_full, binning=binning_next)
+        merged = patched_totals_full.concatenate_bins(appended)
+        binned = patched_totals_full.as_array()[0, 0]
+        npt.assert_array_equal(
+            merged.as_array()[0, 0], np.concatenate([binned, binned])
+        )
+        with raises(ValueError):
+            patched_totals_full.concatenate_bins(patched_totals_full)
+
 
 def patched_counts_from_matrix(binning, matrix, auto):
     n_bins = len(binning)
@@ -307,6 +382,13 @@ def patched_counts_from_matrix(binning, matrix, auto):
     return paircounts.PatchedCount(binning, counts, auto=auto)
 
 
+@fixture
+def patched_counts_full(binning, patch_matrix_full):
+    n_bins = len(binning)
+    matrix = add_zbin_dimension(patch_matrix_full, n_bins)
+    return paircounts.PatchedCount(binning, matrix, auto=False)
+
+
 class TestPatchedCount:
     def test_init(self, binning, patch_matrix_full):
         n_bins = len(binning)
@@ -328,6 +410,27 @@ def test_init(self, binning, patch_matrix_full):
         # just call once
         repr(counts)
 
+    def test_eq(self, patched_counts_full):
+        assert patched_counts_full == patched_counts_full
+        patched_new = paircounts.PatchedCount(
+            patched_counts_full.get_binning(),
+            patched_counts_full.counts + 1,
+            auto=patched_counts_full.auto,
+        )
+        assert patched_counts_full != patched_new
+        patched_new = paircounts.PatchedCount(
+            patched_counts_full.get_binning(),
+            patched_counts_full.counts,
+            auto=not patched_counts_full.auto,
+        )
+        assert patched_counts_full != patched_new
+        assert patched_counts_full != 1
+
+    def test_size(self, patched_counts_full):
+        np = patched_counts_full.n_patches
+        nb = patched_counts_full.n_bins
+        assert patched_counts_full.size == (np * np * nb)
+
     def test_keys_values(self, binning):
         n_bins = len(binning)
         counts = paircounts.PatchedCount.zeros(binning, 2, auto=False)
@@ -404,11 +507,92 @@ def test_jackknife(
         assert result.data[0] == data
         npt.assert_equal(result.samples[:, 0], jack)
 
+        subset = counts.patches[0]
+        assert subset.sample_sum(config).samples.ndim == 2
+
     def test_bootstrap(self, binning, patch_matrix_full):
         counts = patched_counts_from_matrix(binning, patch_matrix_full, auto=False)
         with raises(NotImplementedError):
             counts.sample_sum(ResamplingConfig(method="bootstrap"))
 
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_bins(self, patched_counts_full, items):
+        subset = patched_counts_full.bins[items]
+        npt.assert_array_equal(
+            subset.counts, apply_slice_ndim(patched_counts_full.counts, items, axis=2)
+        )
+        assert (subset.get_binning() == patched_counts_full.get_binning()[items]).all()
+
+    def test_bins_iter(self, patched_counts_full):
+        for i, _bin in enumerate(patched_counts_full.bins):
+            npt.assert_array_equal(
+                _bin.counts, apply_slice_ndim(patched_counts_full.counts, i, axis=2)
+            )
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_patches(self, patched_counts_full, items):
+        subset = patched_counts_full.patches[items]
+        npt.assert_array_equal(
+            subset.counts,
+            apply_slice_ndim(patched_counts_full.counts, items, axis=(0, 1)),
+        )
+        assert (subset.get_binning() == patched_counts_full.get_binning()).all()
+
+    def test_patches_iter(self, patched_counts_full):
+        for i, _bin in enumerate(patched_counts_full.patches):
+            npt.assert_array_equal(
+                _bin.counts,
+                apply_slice_ndim(patched_counts_full.counts, i, axis=(0, 1)),
+            )
+
+    def test_concatenate_patches(self, patched_counts_full, binning_next):
+        merged = patched_counts_full.concatenate_patches(patched_counts_full)
+        diag = np.diag(patched_counts_full.as_array()[:, :, 0])
+        npt.assert_array_equal(
+            np.diag(merged.as_array()[:, :, 0]), np.concatenate([diag, diag])
+        )
+        mod = modify_counts(patched_counts_full, binning=binning_next)
+        with raises(ValueError):
+            patched_counts_full.concatenate_patches(mod)
+
+    def test_concatenate_bins(self, patched_counts_full, binning_next):
+        appended = modify_counts(patched_counts_full, binning=binning_next)
+        merged = patched_counts_full.concatenate_bins(appended)
+        binned = patched_counts_full.as_array()[0, 0]
+        npt.assert_array_equal(
+            merged.as_array()[0, 0], np.concatenate([binned, binned])
+        )
+        with raises(ValueError):
+            patched_counts_full.concatenate_bins(patched_counts_full)
+
+    def test_add(self, patched_counts_full):
+        summed = patched_counts_full + patched_counts_full
+        npt.assert_array_equal(summed.counts, patched_counts_full.counts * 2)
+        with raises(TypeError):
+            patched_counts_full + 1
+        with raises(ValueError):
+            patched_counts_full + patched_counts_full.patches[:-1]
+        with raises(ValueError):
+            patched_counts_full + patched_counts_full.bins[:-1]
+
+    def test_radd(self, patched_counts_full):
+        assert (0 + patched_counts_full) == patched_counts_full
+        assert sum([patched_counts_full]) == patched_counts_full
+        with raises(TypeError):
+            1 + patched_counts_full
+
+    def test_mul(self, patched_counts_full):
+        assert patched_counts_full * 2 == patched_counts_full + patched_counts_full
+        with raises(TypeError):
+            2 * patched_counts_full
+        with raises(TypeError):
+            patched_counts_full * True
+        with raises(TypeError):
+            patched_counts_full * [1, 2]
+
+    def test_sum(self, patched_counts_full):
+        assert patched_counts_full.sum() == np.sum(patched_counts_full.counts)
+
 
 @fixture
 def pair_count_result(patched_totals_full, patch_matrix_full):
@@ -441,6 +625,11 @@ def test_init(self, patched_totals_full, patch_matrix_full):
         # just call once
         repr(res)
 
+    def test_eq(self, pair_count_result):
+        assert pair_count_result == pair_count_result
+        assert pair_count_result != pair_count_result.bins[1:]
+        assert pair_count_result != 1
+
     def test_sample(self, pair_count_result, expect_matrix_full):
         config = ResamplingConfig(method="jackknife")
         result = pair_count_result.sample(config)
@@ -468,3 +657,89 @@ def test_hdf(self, pair_count_result, tmp_hdf5):
             assert binA.left == binB.left
             assert binA.right == binB.right
         assert restored.total.auto == pair_count_result.total.auto
+
+    def test_add(self, pair_count_result):
+        summed = pair_count_result + pair_count_result
+        assert pair_count_result.total == summed.total
+        assert pair_count_result.count * 2 == summed.count
+        with raises(TypeError):
+            pair_count_result + 1
+        # total does not agree
+        new_total = paircounts.PatchedTotal(
+            pair_count_result.get_binning(),
+            pair_count_result.total.totals1 + 1,
+            pair_count_result.total.totals2,
+            auto=pair_count_result.auto,
+        )
+        new = paircounts.NormalisedCounts(pair_count_result.count, new_total)
+        with raises(ValueError, match="total number"):
+            pair_count_result + new
+
+    def test_radd(self, pair_count_result):
+        assert (0 + pair_count_result) == pair_count_result
+        assert sum([pair_count_result]) == pair_count_result
+        with raises(TypeError):
+            1 + pair_count_result
+
+    def test_mul(self, pair_count_result):
+        mul = pair_count_result * 2
+        assert mul.total == pair_count_result.total
+        assert mul.count == pair_count_result.count * 2
+        with raises(TypeError):
+            pair_count_result * True
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_bins(self, pair_count_result, items):
+        pc = pair_count_result
+        assert pc.bins[items] == paircounts.NormalisedCounts(
+            pc.count.bins[items], pc.total.bins[items]
+        )
+
+    @mark.parametrize("items", [1, 2, slice(0, 2), [1, 2]])
+    def test_patches(self, pair_count_result, items):
+        pc = pair_count_result
+        assert pc.patches[items] == paircounts.NormalisedCounts(
+            pc.count.patches[items], pc.total.patches[items]
+        )
+
+    def test_auto(self, pair_count_result):
+        assert pair_count_result.auto == pair_count_result.count.auto
+
+    def test_concatenate_patches(self, pair_count_result):
+        merged = pair_count_result.concatenate_patches(pair_count_result)
+        merged_count = pair_count_result.count.concatenate_patches(
+            pair_count_result.count
+        )
+        merged_total = pair_count_result.total.concatenate_patches(
+            pair_count_result.total
+        )
+        assert merged == paircounts.NormalisedCounts(
+            count=merged_count, total=merged_total
+        )
+        with raises(TypeError, match="concat"):
+            pair_count_result.concatenate_patches(pair_count_result.total)
+
+    def test_concatenate_bins(self, pair_count_result, binning_next):
+        shifted_count = modify_counts(pair_count_result.count, binning=binning_next)
+        shifted_total = modify_totals(pair_count_result.total, binning=binning_next)
+        shifted_pc = paircounts.NormalisedCounts(
+            count=shifted_count, total=shifted_total
+        )
+        merged = pair_count_result.concatenate_bins(shifted_pc)
+        merged_count = pair_count_result.count.concatenate_bins(shifted_count)
+        merged_total = pair_count_result.total.concatenate_bins(shifted_total)
+        assert merged == paircounts.NormalisedCounts(
+            count=merged_count, total=merged_total
+        )
+        with raises(TypeError, match="concat"):
+            pair_count_result.concatenate_bins(pair_count_result.total)
+
+
+def test_pack_results(pair_count_result):
+    pc = pair_count_result
+    counts = dict(a=pc.count)
+    assert paircounts.pack_results(counts, pc.total) == pc
+    counts["b"] = pc.count
+    result = paircounts.pack_results(counts, pc.total)
+    for key in counts:
+        assert result[key] == pc
diff --git a/tests/yaw/test_redshifts.py b/tests/yaw/test_redshifts.py
new file mode 100644
index 00000000..5e826b4c
--- /dev/null
+++ b/tests/yaw/test_redshifts.py
@@ -0,0 +1,53 @@
+from pytest import mark
+
+
+class TestRedshiftData:
+    @mark.skip
+    def test_from_corrdata(self):
+        assert 0
+
+    @mark.skip
+    def test_from_corrfuncs(self):
+        assert 0
+
+    @mark.skip
+    def test_io(self):
+        assert 0
+
+    @mark.skip
+    def test_normalised(self):
+        assert 0
+
+    @mark.skip
+    def test_mean(self):
+        assert 0
+
+    @mark.skip
+    def test_rebin(self):
+        assert 0
+
+    @mark.skip
+    def test_shift(self):
+        assert 0
+
+
+class TestHisttData:
+    @mark.skip
+    def test_io(self):
+        assert 0
+
+    @mark.skip
+    def test_normalised(self):
+        assert 0
+
+    @mark.skip
+    def test_mean(self):
+        assert 0
+
+    @mark.skip
+    def test_rebin(self):
+        assert 0
+
+    @mark.skip
+    def test_shift(self):
+        assert 0