fix: return all-NaN pycnophylactic output for empty valid masks

xrspatial/dasymetric.py

@@ -644,6 +644,8 @@ def _pycnophylactic_numpy(zones_arr, values_dict, nodata_zone,
 
     # pixels that belong to some zone (valid for smoothing)
     valid = ~np.isnan(surface)
+    if not np.any(valid):
+        return surface
 
     for _ in range(max_iterations):
         # Laplacian smoothing: mean of 4-connected neighbours

xrspatial/tests/test_dasymetric.py

@@ -1070,13 +1070,7 @@ def test_three_class_conservation_only(self):
 # ---------------------------------------------------------------------------
 
 class TestPycnophylacticEmptyValid:
-    """pycnophylactic crashes when no pixel is valid for smoothing (#3406).
-
-    disaggregate handles the same inputs gracefully (all-NaN output); these
-    are xfail(strict) until #3406 makes pycnophylactic agree.  When the
-    source fix lands, the tests start XPASSing and strict mode flips them
-    red, prompting removal of the marker.
-    """
+    """pycnophylactic agrees with disaggregate on empty-valid inputs (#3406)."""
 
     def test_disaggregate_all_nan_zones_is_all_nan(self):
         """Reference: disaggregate returns all-NaN, no crash."""
@@ -1085,21 +1079,11 @@ def test_disaggregate_all_nan_zones_is_all_nan(self):
         result = disaggregate(zones, {1: 100.0}, weight)
         assert np.all(np.isnan(result.values))
 
-    @pytest.mark.xfail(
-        reason="#3406: pycnophylactic raises ValueError on empty-valid input",
-        strict=True,
-        raises=ValueError,
-    )
     def test_pycnophylactic_all_nan_zones(self):
         zones = create_test_raster(np.full((3, 3), np.nan), backend='numpy')
         result = pycnophylactic(zones, {1: 100.0})
         assert np.all(np.isnan(result.values))
 
-    @pytest.mark.xfail(
-        reason="#3406: pycnophylactic raises ValueError on empty-valid input",
-        strict=True,
-        raises=ValueError,
-    )
     def test_pycnophylactic_no_matching_zone(self):
         zones = create_test_raster(
             np.array([[1, 1], [2, 2]], dtype=np.float64), backend='numpy')