Tests/catalog methods

seismostats/catalogs/catalog.py

@@ -328,6 +328,7 @@ def estimate_mc(
                                  self.magnitude.max(),
                                  delta_m)
 
+        # TODO change once we have a global estimate_mc
         mc_est = mc_ks(self.magnitude,
                        mcs_test,
                        delta_m,
@@ -348,7 +349,7 @@ def estimate_b(
         weights: list | None = None,
         b_parameter: str = "b_value",
         return_std: bool = False,
-        method: str = "tinti",
+        method: str = "classic",
         return_n: bool = False,
     ) -> float | tuple[float, float] | tuple[float, float, float]:
         """
@@ -394,14 +395,17 @@ def estimate_b(
         if delta_m is None:
             delta_m = self.delta_m
 
+        # filter magnitudes above mc without changing the original dataframe
+        df = self[self.magnitude >= mc]
+
         if method == "positive":
             # dataframe needs 'time' column to be sorted
-            if 'time' not in self.columns:
+            if 'time' not in df.columns:
                 raise ValueError('"time" column needs to be set in order to \
                                  use the b-positive method')
-            mags = self.sort_values("time").magnitude
+            mags = df.sort_values("time").magnitude
         else:
-            mags = self.magnitude
+            mags = df.magnitude
 
         b_estimate = estimate_b(mags,
                                 mc,

seismostats/catalogs/tests/test_catalog.py

@@ -1,12 +1,15 @@
 import os
 import re
 import uuid
+import pytest
 
+import numpy as np
 import pandas as pd
 
 from seismostats.catalogs.catalog import (REQUIRED_COLS_CATALOG, Catalog,
                                           ForecastCatalog)
 from seismostats.utils.binning import bin_to_precision
+from seismostats.analysis.estimate_beta import estimate_b
 
 RAW_DATA = {'name': ['Object 1', 'Object 2', 'Object 3'],
             'magnitude': [10.0, 12.5, 8.2],
@@ -77,19 +80,70 @@ def test_forecast_catalog_strip():
     assert isinstance(dropped, Catalog)
 
 
-def test_catalog_bin():
-    mag_values = [0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]
-    delta_m = 0.1
-
+@pytest.mark.parametrize(
+    "mag_values, delta_m",
+    [
+        (np.array([0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]),
+         0.1),
+        (np.array([0.235, -0.235, 4.499, 5.5, 6, 0.1, 1.6]),
+         0.2),
+        ([0.235, -0.235, 4.499, 5.5, 6, 0.1, 1.6],
+         0.2)
+    ]
+)
+def test_catalog_bin(mag_values: np.ndarray, delta_m: float):
     catalog = Catalog({'magnitude': mag_values})
 
     assert (catalog.bin_magnitudes(
         delta_m)['magnitude'].tolist()
         == bin_to_precision(mag_values, delta_m)).all()
+    assert (catalog.bin_magnitudes()['magnitude'].tolist()
+            == bin_to_precision(mag_values, 0.1)).all()
 
-    catalog.bin_magnitudes(delta_m, inplace=True)
+    return_value = catalog.bin_magnitudes(delta_m, inplace=True)
     assert (catalog['magnitude'].tolist()
             == bin_to_precision(mag_values, delta_m)).all()
+    assert return_value is None
+    assert catalog.delta_m == delta_m
+
+
+def test_catalog_estimate_mc():
+    # TODO once global_mc method is implemented
+    catalog = Catalog({'magnitude': [0.235, -0.235, 4.499, 4.5, 6, 0.1, 1.6]})
+
+    with pytest.raises(ValueError):
+        catalog.estimate_mc()
+
+
+@pytest.mark.parametrize(
+    "mag_values, delta_m, mc",
+    [
+        (np.array([0.0, 0.235, 0.238, 4.499, 4.5, 6, 0.1, 1.6]),
+         0.001, 0.0)
+    ]
+)
+def test_catalog_estimate_b(mag_values, delta_m, mc):
+    catalog = Catalog({'magnitude': mag_values})
+
+    with pytest.raises(ValueError):
+        catalog.estimate_b(mc=None, delta_m=None)
+        catalog.estimate_b(mc=1.0, delta_m=None)
+        catalog.estimate_b(mc=None, delta_m=0.1)
+        catalog.estimate_b(mc=1.0, delta_m=0.1, method='positive')
+
+    b_value = estimate_b(catalog['magnitude'],
+                         mc=mc,
+                         delta_m=delta_m,
+                         method="classic")
+    return_value = catalog.estimate_b(mc=mc, delta_m=delta_m, method="classic")
+    assert catalog.b_value == b_value
+    assert return_value == b_value
+
+    catalog.mc = mc
+    catalog.delta_m = delta_m
+    return_value = catalog.estimate_b(method="classic")
+    assert catalog.b_value == b_value
+    assert return_value == b_value
 
 
 def test_to_quakeml():

seismostats/utils/binning.py

@@ -52,9 +52,12 @@ def bin_to_precision(x: np.ndarray | list, delta_x: float = 0.1) -> np.ndarray:
     """
     if x is None:
         raise ValueError("x cannot be None")
+    if delta_x == 0:
+        raise ValueError("delta_x cannot be 0")
 
     if isinstance(x, list):
         x = np.array(x)
+
     d = decimal.Decimal(str(delta_x))
     decimal_places = abs(d.as_tuple().exponent)
     return np.round(normal_round_to_int(x / delta_x) * delta_x, decimal_places)

seismostats/utils/tests/test_binning.py

@@ -43,10 +43,15 @@ def test_bin_to_precision(x: np.ndarray, delta_x: float,
     y = bin_to_precision(x, delta_x)
     assert (y == rounded_value).all()
 
+    y = bin_to_precision(x)
+    z = bin_to_precision(x, 0.1)
+    assert (y == z).all()
+
 
 def test_bin_to_precision_none():
     with pytest.raises(ValueError):
         bin_to_precision(None)
+        bin_to_precision([], 0)
 
 
 @pytest.mark.parametrize(