diff --git a/docs/source/tutorials/gaussian_laser.ipynb b/docs/source/tutorials/gaussian_laser.ipynb
index 850bfa26..e9b99f20 100644
--- a/docs/source/tutorials/gaussian_laser.ipynb
+++ b/docs/source/tutorials/gaussian_laser.ipynb
@@ -78,6 +78,24 @@
     "laser = Laser(dimensions,lo,hi,num_points,laser_profile)"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "8ac8f985-4a09-4b9a-a9ff-11ce2ef94caa",
+   "metadata": {},
+   "source": [
+    "The laser pulse can be visualized with the `show` method."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "bb8dc634-b3c4-4448-834b-fe9acb6c9645",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "laser.show()"
+   ]
+  },
   {
    "cell_type": "markdown",
    "id": "292c6e63-0480-4fb9-b1ad-6b679a3472d0",
@@ -97,6 +115,16 @@
     "laser.propagate(-z_R)                               # Propagate the pulse upstream of the focal plane"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "3158e51d-e331-438d-90c1-5f8c63bf9841",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "laser.show()"
+   ]
+  },
   {
    "cell_type": "markdown",
    "id": "8df816dc-9b47-4cb0-8716-600352fafb72",
@@ -143,7 +171,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.5"
+   "version": "3.11.6"
   }
  },
  "nbformat": 4,
diff --git a/lasy/laser.py b/lasy/laser.py
index 7b7503ea..f7eff8a5 100644
--- a/lasy/laser.py
+++ b/lasy/laser.py
@@ -276,3 +276,46 @@ def write_to_file(
             self.profile.pol,
             save_as_vector_potential,
         )
+
+    def show(self, **kw):
+        """
+        Show a 2D image of the laser amplitude.
+
+        Parameters
+        ----------
+        **kw: additional arguments to be passed to matplotlib's imshow command
+        """
+        if self.dim == "rt":
+            # Show field in the plane y=0, above and below axis, with proper sign for each mode
+            E = [
+                np.concatenate(
+                    ((-1) ** m * self.grid.field[0, ::-1], self.grid.field[0])
+                )
+                for m in self.grid.azimuthal_modes
+            ]
+            E = sum(E)  # Sum all the modes
+            extent = [
+                self.grid.lo[-1],
+                self.grid.hi[-1],
+                -self.grid.hi[0],
+                self.grid.hi[0],
+            ]
+
+        else:
+            # In 3D show an image in the xt plane
+            i_slice = int(self.grid.field.shape[1] // 2)
+            E = self.grid.field[:, i_slice, :]
+            extent = [
+                self.grid.lo[-1],
+                self.grid.hi[-1],
+                self.grid.lo[0],
+                self.grid.hi[0],
+            ]
+
+        import matplotlib.pyplot as plt
+
+        plt.imshow(abs(E), extent=extent, aspect="auto", origin="lower", **kw)
+        cb = plt.colorbar()
+        cb.set_label("$|E_{envelope}|$ (V/m)")
+        plt.xlabel("t (s)")
+        plt.ylabel("x (m)")