diff --git a/python/simulation.py b/python/simulation.py
index 42e257a61..9afcdb94c 100644
--- a/python/simulation.py
+++ b/python/simulation.py
@@ -4918,7 +4918,7 @@ def plot_fields(self, **kwargs):
         return vis.plot_fields(self, **kwargs)
 
     def plot3D(
-        self, save_to_image: bool = False, image_name: str = "sim.png", **kwargs
+        self, save_to_image: bool = False, image_name: str = "sim.png", eps_frequency: float = 0., resolution: float | None = None, **kwargs
     ):
         """
         Uses vispy to render a 3D scene of the simulation object. The simulation object must be 3D.
@@ -4932,10 +4932,16 @@ def plot3D(
             scale_factor: float, camera zoom factor
             azimuth: float, azimuthal angle in degrees
             elevation: float, elevation angle in degrees
+            eps_frequency: for materials with a [frequency-dependent
+                        permittivity](Materials.md#material-dispersion) $\\varepsilon(f)$, specifies the
+                        frequency $f$ (in Meep units) of the real part of the permittivity to use in the
+                        plot. Defaults to 0.
+            resolution: the resolution used for plotting. Defaults to the
+                        `resolution` of the `Simulation` object.
         """
         import meep.visualization as vis
 
-        return vis.plot3D(self, save_to_image, image_name, **kwargs)
+        return vis.plot3D(self, save_to_image, image_name, eps_frequency, resolution, **kwargs)
 
     def visualize_chunks(self):
         """
diff --git a/python/visualization.py b/python/visualization.py
index 55c8a04e8..df83af882 100644
--- a/python/visualization.py
+++ b/python/visualization.py
@@ -1089,7 +1089,7 @@ def plot2D(
 CELL_EDGE_COLOR_3D: tuple[float, float, float, float] = (0.75, 0.75, 0.75, 1)  # gray
 
 
-def plot3D(sim, save_to_image: bool = False, image_name: str = "sim.png", **kwargs):
+def plot3D(sim: mp.Simulation, save_to_image: bool = False, image_name: str = "sim.png", eps_frequency: float = 0, resolution: float | None = None, **kwargs):
     from vispy.scene.visuals import Box, Mesh
     from vispy.scene import SceneCanvas, transforms
 
@@ -1112,7 +1112,8 @@ def plot3D(sim, save_to_image: bool = False, image_name: str = "sim.png", **kwar
         sim_center, sim_size, sim.is_cylindrical
     )
 
-    grid_resolution = sim.resolution
+    # Get eps parameters or use default
+    grid_resolution = resolution or sim.resolution
 
     Nx = int((xmax - xmin) * grid_resolution + 1)
     Ny = int((ymax - ymin) * grid_resolution + 1)
@@ -1123,12 +1124,13 @@ def plot3D(sim, save_to_image: bool = False, image_name: str = "sim.png", **kwar
     ztics = np.linspace(zmin, zmax, Nz)
 
     # Get eps for geometry
-    eps_data = np.round(np.real(sim.get_epsilon_grid(xtics, ytics, ztics)), 2)
 
-    unique = np.unique(np.abs(eps_data)).tolist()
+    eps_data = np.real(sim.get_epsilon_grid(xtics, ytics, ztics, eps_frequency))
+
+    unique = np.unique(eps_data).tolist()
 
     # Remove background material
-    unique.remove(np.round(np.abs(np.asarray(sim.default_material.epsilon_diag)), 2)[0])
+    unique.remove(np.asarray(sim.default_material.epsilon_diag)[0])
 
     mesh_midpoint = (sim_size[0] / 2, sim_size[1] / 2, sim_size[2] / 2)