Fix NaN values in inverse Box-Cox transform for Lambda > 0

pysteps/tests/test_utils_transformation.py

@@ -84,6 +84,24 @@ def test_boxcox_transform(R, metadata, Lambda, threshold, zerovalue, inverse, ex
     )
 
 
+@pytest.mark.parametrize("Lambda", [0.0, 0.5, 1.0])
+def test_boxcox_transform_roundtrip_with_zeros(Lambda):
+    """Test that the inverse boxcox_transform maps the zero value back to zero."""
+    R = np.array([0.0, 0.05, 0.1, 1.0, 10.0])
+    metadata = {
+        "accutime": 5,
+        "transform": None,
+        "unit": "mm/h",
+        "threshold": 0.1,
+        "zerovalue": 0.0,
+    }
+    R_t, metadata_t = transformation.boxcox_transform(R.copy(), metadata, Lambda)
+    R_b, metadata_b = transformation.boxcox_transform(R_t, metadata_t, inverse=True)
+    assert_array_almost_equal(R_b, [0.0, 0.0, 0.1, 1.0, 10.0])
+    assert metadata_b["transform"] is None
+    assert metadata_b["threshold"] == pytest.approx(0.1)
+
+
 # dB_transform
 test_data = [
     (

pysteps/utils/transformation.py

@@ -129,6 +129,10 @@ def boxcox_transform(
         if zerovalue is None:
             zerovalue = 0.0
 
+        # Find the zeros before the back-transformation: for Lambda > 0, values
+        # below -1 / Lambda (e.g. the zero value) have no inverse and become NaN
+        zeros = R < threshold
+
         # Apply inverse Box-Cox transform
         if Lambda == 0.0:
             R = np.exp(R)
@@ -138,7 +142,7 @@ def boxcox_transform(
             R = np.exp(np.log(Lambda * R + 1) / Lambda)
             threshold = np.exp(np.log(Lambda * threshold + 1) / Lambda)
 
-        R[R < threshold] = zerovalue
+        R[zeros] = zerovalue
 
         metadata["transform"] = None
         metadata["zerovalue"] = zerovalue