diff --git a/xrspatial/rasterize.py b/xrspatial/rasterize.py
index 71e74f68f..0dc381e90 100644
--- a/xrspatial/rasterize.py
+++ b/xrspatial/rasterize.py
@@ -1951,6 +1951,11 @@ def _parse_input(geometries, column=None, columns=None):
         if isinstance(geometries, gpd.GeoDataFrame):
             geom_list = geometries.geometry.tolist()
             total_bounds = tuple(geometries.total_bounds)
+            # GeoPandas returns (nan, nan, nan, nan) for an empty frame.
+            # Treat as "no inferred bounds" so the caller's explicit-bounds
+            # guard fires instead of producing a raster with nan coords.
+            if any(not np.isfinite(v) for v in total_bounds):
+                total_bounds = None
             if columns is not None:
                 props_array = geometries[columns].values.astype(np.float64)
             else:
@@ -2217,19 +2222,28 @@ def rasterize(
     if final_bounds is None:
         final_bounds = inferred_bounds
     if final_bounds is None and geom_list:
-        # Compute bounds lazily only when not supplied by the caller
+        # Compute bounds lazily only when not supplied by the caller.
+        # Drop empty/invalid geoms whose bbox is nan so they don't poison
+        # the inferred extent.
         geom_bboxes = _geometry_bboxes(geom_list)
         if len(geom_bboxes) > 0:
-            final_bounds = (geom_bboxes[:, 0].min(),
-                            geom_bboxes[:, 1].min(),
-                            geom_bboxes[:, 2].max(),
-                            geom_bboxes[:, 3].max())
+            finite = np.all(np.isfinite(geom_bboxes), axis=1)
+            geom_bboxes = geom_bboxes[finite]
+            if len(geom_bboxes) > 0:
+                final_bounds = (geom_bboxes[:, 0].min(),
+                                geom_bboxes[:, 1].min(),
+                                geom_bboxes[:, 2].max(),
+                                geom_bboxes[:, 3].max())
     if final_bounds is None:
         raise ValueError(
             "bounds must be provided when geometries are empty or have "
             "no spatial extent")
 
     xmin, ymin, xmax, ymax = final_bounds
+    if not all(np.isfinite(v) for v in (xmin, ymin, xmax, ymax)):
+        raise ValueError(
+            f"Invalid bounds: all of (xmin, ymin, xmax, ymax) must be "
+            f"finite, got {(xmin, ymin, xmax, ymax)!r}")
     if xmin >= xmax or ymin >= ymax:
         raise ValueError(
             f"Invalid bounds: xmin ({xmin}) must be < xmax ({xmax}) and "
diff --git a/xrspatial/tests/test_rasterize.py b/xrspatial/tests/test_rasterize.py
index 20d3a33b8..78a993ba7 100644
--- a/xrspatial/tests/test_rasterize.py
+++ b/xrspatial/tests/test_rasterize.py
@@ -364,6 +364,49 @@ def test_invalid_chunks(self, bad):
                       width=10, height=10, bounds=(0, 0, 1, 1),
                       chunks=bad)
 
+    def test_empty_geometry_does_not_poison_inferred_bounds(self):
+        # Issue #2065: an empty geometry returns nan from .bounds; the
+        # caller used .min()/.max() unfiltered, so a single empty geom
+        # poisoned the inferred extent and produced a raster with nan
+        # x/y coords.  Drop the empty and infer from the rest.
+        empty = Polygon()
+        result = rasterize(
+            [(empty, 99), (box(0, 0, 1, 1), 1)],
+            width=2, height=2)
+        assert np.all(np.isfinite(result.x.values))
+        assert np.all(np.isfinite(result.y.values))
+        assert np.all(result.values == 1)
+
+    def test_only_empty_geometry_requires_explicit_bounds(self):
+        empty = Polygon()
+        with pytest.raises(ValueError, match="bounds must be provided"):
+            rasterize([(empty, 1.0)], width=2, height=2)
+
+    def test_empty_geodataframe_requires_explicit_bounds(self):
+        # Issue #2065: total_bounds on an empty frame is (nan,nan,nan,nan).
+        # That used to bypass the empty-bounds guard.
+        import geopandas as gpd
+        empty_gdf = gpd.GeoDataFrame(
+            {'value': []}, geometry=gpd.GeoSeries([]))
+        with pytest.raises(ValueError, match="bounds must be provided"):
+            rasterize(empty_gdf, width=2, height=2, column='value')
+
+    def test_explicit_nan_bounds_rejected(self):
+        with pytest.raises(ValueError, match="must be finite"):
+            rasterize([(box(0, 0, 1, 1), 1.0)], width=2, height=2,
+                      bounds=(0, 0, float('nan'), 1))
+
+    def test_empty_multipolygon_filtered_from_inferred_bounds(self):
+        # Same shape as the empty-Polygon case but exercises the
+        # MultiPolygon branch in shapely's is_empty / bounds path.
+        empty_mp = MultiPolygon()
+        result = rasterize(
+            [(empty_mp, 99), (box(0, 0, 1, 1), 1)],
+            width=2, height=2)
+        assert np.all(np.isfinite(result.x.values))
+        assert np.all(np.isfinite(result.y.values))
+        assert np.all(result.values == 1)
+
 
 # ---------------------------------------------------------------------------
 # all_touched mode