diff --git a/autolens/analysis/plotter.py b/autolens/analysis/plotter.py
index 67a07f032..7a31868a9 100644
--- a/autolens/analysis/plotter.py
+++ b/autolens/analysis/plotter.py
@@ -72,6 +72,7 @@ def should_plot(name):
                 grid=grid,
                 output_path=output_path,
                 output_format=fmt,
+                title_prefix=self.title_prefix,
             )
 
         if should_plot("fits_tracer"):
@@ -104,10 +105,11 @@ def should_plot(name):
             if isinstance(fmt, (list, tuple)):
                 fmt = fmt[0]
 
+            _title = f"{self.title_prefix.rstrip()} Image With Positions" if self.title_prefix else "Image With Positions"
             plot_array(
                 array=image,
                 positions=[pos_arr],
-                title="Image With Positions",
+                title=_title,
                 output_path=str(self.image_path),
                 output_filename="image_with_positions",
                 output_format=fmt,
diff --git a/autolens/imaging/model/plotter.py b/autolens/imaging/model/plotter.py
index 5fdeb1b92..061cbad73 100644
--- a/autolens/imaging/model/plotter.py
+++ b/autolens/imaging/model/plotter.py
@@ -84,12 +84,14 @@ def should_plot(name):
                         plane_index=plane_index,
                         image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                         source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                        title_prefix=self.title_prefix,
                     )
             else:
                 subplot_fit(
                     fit, output_path=output_path, output_format=fmt,
                     image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                     source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                    title_prefix=self.title_prefix,
                 )
 
         if quick_update:
@@ -100,6 +102,7 @@ def should_plot(name):
                 fit, output_path=output_path, output_format=fmt,
                 image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                 source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                title_prefix=self.title_prefix,
             )
 
         if should_plot("subplot_fit_log10"):
@@ -111,18 +114,21 @@ def should_plot(name):
                             plane_index=plane_index,
                             image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                             source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                            title_prefix=self.title_prefix,
                         )
                 else:
                     subplot_fit_log10(
                         fit, output_path=output_path, output_format=fmt,
                         image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                         source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                        title_prefix=self.title_prefix,
                     )
             except ValueError:
                 pass
 
         if should_plot("subplot_of_planes"):
-            subplot_of_planes(fit, output_path=output_path, output_format=fmt)
+            subplot_of_planes(fit, output_path=output_path, output_format=fmt,
+                              title_prefix=self.title_prefix)
 
         if should_plot("fits_fit"):
             fits_fit(fit=fit, output_path=self.image_path)
@@ -154,9 +160,11 @@ def should_plot(name):
         fmt = self.fmt
 
         if should_plot("subplot_fit") or quick_update:
-            subplot_fit_combined(fit_list, output_path=output_path, output_format=fmt)
+            subplot_fit_combined(fit_list, output_path=output_path, output_format=fmt,
+                                 title_prefix=self.title_prefix)
 
             if quick_update:
                 return
 
-            subplot_fit_combined_log10(fit_list, output_path=output_path, output_format=fmt)
+            subplot_fit_combined_log10(fit_list, output_path=output_path, output_format=fmt,
+                                       title_prefix=self.title_prefix)
diff --git a/autolens/imaging/plot/fit_imaging_plots.py b/autolens/imaging/plot/fit_imaging_plots.py
index d6598964c..a9f7ce069 100644
--- a/autolens/imaging/plot/fit_imaging_plots.py
+++ b/autolens/imaging/plot/fit_imaging_plots.py
@@ -189,6 +189,7 @@ def subplot_fit(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 12-panel subplot summarising an imaging fit.
@@ -228,7 +229,8 @@ def subplot_fit(
     """
     if len(fit.tracer.planes) == 1:
         return subplot_fit_x1_plane(fit, output_path=output_path,
-                                    output_format=output_format, colormap=colormap)
+                                    output_format=output_format, colormap=colormap,
+                                    title_prefix=title_prefix)
 
     plane_index_tag = "" if plane_index is None else f"_{plane_index}"
     final_plane_index = (
@@ -242,18 +244,19 @@ def subplot_fit(
             _compute_critical_curves_from_fit(fit)
         )
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 4, figsize=conf_subplot_figsize(3, 4))
     axes_flat = list(axes.flatten())
 
-    plot_array(array=fit.data, ax=axes_flat[0], title="Data", colormap=colormap)
+    plot_array(array=fit.data, ax=axes_flat[0], title=_pf("Data"), colormap=colormap)
 
     # Data at source scale
-    plot_array(array=fit.data, ax=axes_flat[1], title="Data (Source Scale)",
+    plot_array(array=fit.data, ax=axes_flat[1], title=_pf("Data (Source Scale)"),
                colormap=colormap, vmax=source_vmax)
 
     plot_array(array=fit.signal_to_noise_map, ax=axes_flat[2],
-               title="Signal-To-Noise Map", colormap=colormap)
-    plot_array(array=fit.model_data, ax=axes_flat[3], title="Model Image",
+               title=_pf("Signal-To-Noise Map"), colormap=colormap)
+    plot_array(array=fit.model_data, ax=axes_flat[3], title=_pf("Model Image"),
                colormap=colormap, lines=image_plane_lines,
                line_colors=image_plane_line_colors)
 
@@ -264,7 +267,7 @@ def subplot_fit(
         lens_model_img = None
     if lens_model_img is not None:
         plot_array(array=lens_model_img, ax=axes_flat[4],
-                   title="Lens Light Model Image", colormap=colormap)
+                   title=_pf("Lens Light Model Image"), colormap=colormap)
     else:
         axes_flat[4].axis("off")
 
@@ -274,7 +277,7 @@ def subplot_fit(
     except (IndexError, AttributeError):
         subtracted_img = None
     if subtracted_img is not None:
-        plot_array(array=subtracted_img, ax=axes_flat[5], title="Lens Light Subtracted",
+        plot_array(array=subtracted_img, ax=axes_flat[5], title=_pf("Lens Light Subtracted"),
                    colormap=colormap, vmin=0.0 if source_vmax is not None else None,
                    vmax=source_vmax)
     else:
@@ -286,7 +289,7 @@ def subplot_fit(
     except (IndexError, AttributeError):
         source_model_img = None
     if source_model_img is not None:
-        plot_array(array=source_model_img, ax=axes_flat[6], title="Source Model Image",
+        plot_array(array=source_model_img, ax=axes_flat[6], title=_pf("Source Model Image"),
                    colormap=colormap, vmax=source_vmax, lines=image_plane_lines,
                    line_colors=image_plane_line_colors)
     else:
@@ -294,27 +297,27 @@ def subplot_fit(
 
     # Source plane zoomed
     _plot_source_plane(fit, axes_flat[7], final_plane_index, zoom_to_brightest=True,
-                       colormap=colormap, title="Source Plane (Zoomed)",
+                       colormap=colormap, title=_pf("Source Plane (Zoomed)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors,
                        vmax=source_vmax)
 
     # Normalized residual map (symmetric)
     norm_resid = fit.normalized_residual_map
     _abs_max = _symmetric_vmax(norm_resid)
-    plot_array(array=norm_resid, ax=axes_flat[8], title="Normalized Residual Map",
+    plot_array(array=norm_resid, ax=axes_flat[8], title=_pf("Normalized Residual Map"),
                colormap=colormap, vmin=-_abs_max, vmax=_abs_max)
 
     # Normalized residual map clipped to [-1, 1]
     plot_array(array=norm_resid, ax=axes_flat[9],
-               title=r"Normalized Residual Map $1\sigma$",
+               title=_pf(r"Normalized Residual Map $1\sigma$"),
                colormap=colormap, vmin=-1.0, vmax=1.0)
 
     plot_array(array=fit.chi_squared_map, ax=axes_flat[10],
-               title="Chi-Squared Map", colormap=colormap, cb_unit=r"$\chi^2$")
+               title=_pf("Chi-Squared Map"), colormap=colormap, cb_unit=r"$\chi^2$")
 
     # Source plane not zoomed
     _plot_source_plane(fit, axes_flat[11], final_plane_index, zoom_to_brightest=False,
-                       colormap=colormap, title="Source Plane (No Zoom)",
+                       colormap=colormap, title=_pf("Source Plane (No Zoom)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors,
                        vmax=source_vmax)
 
@@ -328,6 +331,7 @@ def subplot_fit_x1_plane(
     output_path: Optional[str] = None,
     output_format: str = None,
     colormap: Optional[str] = None,
+    title_prefix: str = None,
 ):
     """
     Produce a 6-panel subplot for a single-plane tracer imaging fit.
@@ -356,6 +360,7 @@ def subplot_fit_x1_plane(
     colormap : str, optional
         Matplotlib colormap name applied to all image panels.
     """
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(2, 3, figsize=conf_subplot_figsize(2, 3))
     axes_flat = list(axes.flatten())
 
@@ -364,23 +369,23 @@ def subplot_fit_x1_plane(
     except (IndexError, AttributeError, ValueError):
         vmax = None
 
-    plot_array(array=fit.data, ax=axes_flat[0], title="Data", colormap=colormap, vmax=vmax)
+    plot_array(array=fit.data, ax=axes_flat[0], title=_pf("Data"), colormap=colormap, vmax=vmax)
 
     plot_array(array=fit.signal_to_noise_map, ax=axes_flat[1],
-               title="Signal-To-Noise Map", colormap=colormap)
+               title=_pf("Signal-To-Noise Map"), colormap=colormap)
 
-    plot_array(array=fit.model_data, ax=axes_flat[2], title="Model Image",
+    plot_array(array=fit.model_data, ax=axes_flat[2], title=_pf("Model Image"),
                colormap=colormap, vmax=vmax)
 
     norm_resid = fit.normalized_residual_map
-    plot_array(array=norm_resid, ax=axes_flat[3], title="Lens Light Subtracted",
+    plot_array(array=norm_resid, ax=axes_flat[3], title=_pf("Lens Light Subtracted"),
                colormap=colormap, cb_unit=r"$\sigma$")
 
-    plot_array(array=norm_resid, ax=axes_flat[4], title="Subtracted Image Zero Minimum",
+    plot_array(array=norm_resid, ax=axes_flat[4], title=_pf("Subtracted Image Zero Minimum"),
                colormap=colormap, vmin=0.0, cb_unit=r"$\sigma$")
 
     _abs_max = _symmetric_vmax(norm_resid)
-    plot_array(array=norm_resid, ax=axes_flat[5], title="Normalized Residual Map",
+    plot_array(array=norm_resid, ax=axes_flat[5], title=_pf("Normalized Residual Map"),
                colormap=colormap, vmin=-_abs_max, vmax=_abs_max, cb_unit=r"$\sigma$")
 
     tight_layout()
@@ -397,6 +402,7 @@ def subplot_fit_log10(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 12-panel subplot summarising an imaging fit with log10 colour scaling.
@@ -427,7 +433,8 @@ def subplot_fit_log10(
     """
     if len(fit.tracer.planes) == 1:
         return subplot_fit_log10_x1_plane(fit, output_path=output_path,
-                                          output_format=output_format, colormap=colormap)
+                                          output_format=output_format, colormap=colormap,
+                                          title_prefix=title_prefix)
 
     plane_index_tag = "" if plane_index is None else f"_{plane_index}"
     final_plane_index = (
@@ -441,69 +448,72 @@ def subplot_fit_log10(
             _compute_critical_curves_from_fit(fit)
         )
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 4, figsize=conf_subplot_figsize(3, 4))
     axes_flat = list(axes.flatten())
 
-    plot_array(array=fit.data, ax=axes_flat[0], title="Data", colormap=colormap,
+    plot_array(array=fit.data, ax=axes_flat[0], title=_pf("Data"), colormap=colormap,
                use_log10=True)
 
     try:
-        plot_array(array=fit.data, ax=axes_flat[1], title="Data (Source Scale)",
+        plot_array(array=fit.data, ax=axes_flat[1], title=_pf("Data (Source Scale)"),
                    colormap=colormap, use_log10=True)
     except ValueError:
         axes_flat[1].axis("off")
 
     try:
         plot_array(array=fit.signal_to_noise_map, ax=axes_flat[2],
-                   title="Signal-To-Noise Map", colormap=colormap, use_log10=True)
+                   title=_pf("Signal-To-Noise Map"), colormap=colormap, use_log10=True)
     except ValueError:
         axes_flat[2].axis("off")
 
-    plot_array(array=fit.model_data, ax=axes_flat[3], title="Model Image",
+    plot_array(array=fit.model_data, ax=axes_flat[3], title=_pf("Model Image"),
                colormap=colormap, use_log10=True, lines=image_plane_lines,
                line_colors=image_plane_line_colors)
 
     try:
         lens_model_img = fit.model_images_of_planes_list[0]
         plot_array(array=lens_model_img, ax=axes_flat[4],
-                   title="Lens Light Model Image", colormap=colormap, use_log10=True)
+                   title=_pf("Lens Light Model Image"), colormap=colormap, use_log10=True)
     except (IndexError, AttributeError):
         axes_flat[4].axis("off")
 
     try:
         subtracted_img = fit.subtracted_images_of_planes_list[final_plane_index]
         plot_array(array=subtracted_img, ax=axes_flat[5],
-                   title="Lens Light Subtracted", colormap=colormap, use_log10=True)
+                   title=_pf("Lens Light Subtracted"), colormap=colormap, use_log10=True)
     except (IndexError, AttributeError):
         axes_flat[5].axis("off")
 
     try:
         source_model_img = fit.model_images_of_planes_list[final_plane_index]
         plot_array(array=source_model_img, ax=axes_flat[6],
-                   title="Source Model Image", colormap=colormap, use_log10=True,
+                   title=_pf("Source Model Image"), colormap=colormap, use_log10=True,
                    lines=image_plane_lines, line_colors=image_plane_line_colors)
     except (IndexError, AttributeError):
         axes_flat[6].axis("off")
 
     _plot_source_plane(fit, axes_flat[7], final_plane_index, zoom_to_brightest=True,
                        colormap=colormap, use_log10=True,
+                       title=_pf("Source Plane (Zoomed)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors,
                        vmax=source_vmax)
 
     norm_resid = fit.normalized_residual_map
     _abs_max = _symmetric_vmax(norm_resid)
-    plot_array(array=norm_resid, ax=axes_flat[8], title="Normalized Residual Map",
+    plot_array(array=norm_resid, ax=axes_flat[8], title=_pf("Normalized Residual Map"),
                colormap=colormap, vmin=-_abs_max, vmax=_abs_max, cb_unit=r"$\sigma$")
 
     plot_array(array=norm_resid, ax=axes_flat[9],
-               title=r"Normalized Residual Map $1\sigma$",
+               title=_pf(r"Normalized Residual Map $1\sigma$"),
                colormap=colormap, vmin=-1.0, vmax=1.0, cb_unit=r"$\sigma$")
 
-    plot_array(array=fit.chi_squared_map, ax=axes_flat[10], title="Chi-Squared Map",
+    plot_array(array=fit.chi_squared_map, ax=axes_flat[10], title=_pf("Chi-Squared Map"),
                colormap=colormap, use_log10=True, cb_unit=r"$\chi^2$")
 
     _plot_source_plane(fit, axes_flat[11], final_plane_index, zoom_to_brightest=False,
                        colormap=colormap, use_log10=True,
+                       title=_pf("Source Plane (No Zoom)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors,
                        vmax=source_vmax)
 
@@ -516,6 +526,7 @@ def subplot_fit_log10_x1_plane(
     output_path: Optional[str] = None,
     output_format: str = None,
     colormap: Optional[str] = None,
+    title_prefix: str = None,
 ):
     """
     Produce a 6-panel log10 subplot for a single-plane tracer imaging fit.
@@ -539,6 +550,7 @@ def subplot_fit_log10_x1_plane(
     colormap : str, optional
         Matplotlib colormap name applied to all image panels.
     """
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(2, 3, figsize=conf_subplot_figsize(2, 3))
     axes_flat = list(axes.flatten())
 
@@ -547,25 +559,25 @@ def subplot_fit_log10_x1_plane(
     except (IndexError, AttributeError, ValueError):
         vmax = None
 
-    plot_array(array=fit.data, ax=axes_flat[0], title="Data", colormap=colormap,
+    plot_array(array=fit.data, ax=axes_flat[0], title=_pf("Data"), colormap=colormap,
                vmax=vmax, use_log10=True)
 
     try:
         plot_array(array=fit.signal_to_noise_map, ax=axes_flat[1],
-                   title="Signal-To-Noise Map", colormap=colormap, use_log10=True)
+                   title=_pf("Signal-To-Noise Map"), colormap=colormap, use_log10=True)
     except ValueError:
         axes_flat[1].axis("off")
 
-    plot_array(array=fit.model_data, ax=axes_flat[2], title="Model Image",
+    plot_array(array=fit.model_data, ax=axes_flat[2], title=_pf("Model Image"),
                colormap=colormap, vmax=vmax, use_log10=True)
 
     norm_resid = fit.normalized_residual_map
-    plot_array(array=norm_resid, ax=axes_flat[3], title="Lens Light Subtracted",
+    plot_array(array=norm_resid, ax=axes_flat[3], title=_pf("Lens Light Subtracted"),
                colormap=colormap, cb_unit=r"$\sigma$")
     _abs_max = _symmetric_vmax(norm_resid)
-    plot_array(array=norm_resid, ax=axes_flat[4], title="Normalized Residual Map",
+    plot_array(array=norm_resid, ax=axes_flat[4], title=_pf("Normalized Residual Map"),
                colormap=colormap, vmin=-_abs_max, vmax=_abs_max, cb_unit=r"$\sigma$")
-    plot_array(array=fit.chi_squared_map, ax=axes_flat[5], title="Chi-Squared Map",
+    plot_array(array=fit.chi_squared_map, ax=axes_flat[5], title=_pf("Chi-Squared Map"),
                colormap=colormap, use_log10=True, cb_unit=r"$\chi^2$")
 
     tight_layout()
@@ -578,6 +590,7 @@ def subplot_of_planes(
     output_format: str = None,
     colormap: Optional[str] = None,
     plane_index: Optional[int] = None,
+    title_prefix: str = None,
 ):
     """
     Produce a 4-panel subplot for each plane in the tracer.
@@ -614,27 +627,29 @@ def subplot_of_planes(
     else:
         plane_indexes = [plane_index]
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     for pidx in plane_indexes:
         fig, axes = subplots(1, 4, figsize=conf_subplot_figsize(1, 4))
         axes_flat = list(axes.flatten())
 
-        plot_array(array=fit.data, ax=axes_flat[0], title="Data", colormap=colormap)
+        plot_array(array=fit.data, ax=axes_flat[0], title=_pf("Data"), colormap=colormap)
 
         try:
             subtracted = fit.subtracted_images_of_planes_list[pidx]
             plot_array(array=subtracted, ax=axes_flat[1],
-                       title=f"Subtracted Image Plane {pidx}", colormap=colormap)
+                       title=_pf(f"Subtracted Image Plane {pidx}"), colormap=colormap)
         except (IndexError, AttributeError):
             axes_flat[1].axis("off")
 
         try:
             model_img = fit.model_images_of_planes_list[pidx]
             plot_array(array=model_img, ax=axes_flat[2],
-                       title=f"Model Image Plane {pidx}", colormap=colormap)
+                       title=_pf(f"Model Image Plane {pidx}"), colormap=colormap)
         except (IndexError, AttributeError):
             axes_flat[2].axis("off")
 
-        _plot_source_plane(fit, axes_flat[3], pidx, colormap=colormap)
+        _plot_source_plane(fit, axes_flat[3], pidx, colormap=colormap,
+                           title=_pf(f"Source Plane {pidx}"))
 
         tight_layout()
         save_figure(fig, path=output_path, filename=f"fit_of_plane_{pidx}", format=output_format)
@@ -649,6 +664,7 @@ def subplot_tracer_from_fit(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 9-panel tracer subplot derived from a `FitImaging` object.
@@ -701,11 +717,12 @@ def subplot_tracer_from_fit(
 
     magnification = LensCalc.from_mass_obj(tracer).magnification_2d_from(grid=grid)
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 3, figsize=conf_subplot_figsize(3, 3))
     axes_flat = list(axes.flatten())
 
     # Panel 0: Model Image
-    plot_array(array=fit.model_data, ax=axes_flat[0], title="Model Image",
+    plot_array(array=fit.model_data, ax=axes_flat[0], title=_pf("Model Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                colormap=colormap)
 
@@ -715,28 +732,27 @@ def subplot_tracer_from_fit(
     except (IndexError, AttributeError):
         source_model_img = None
     if source_model_img is not None:
-        plot_array(array=source_model_img, ax=axes_flat[1], title="Source Model Image",
+        plot_array(array=source_model_img, ax=axes_flat[1], title=_pf("Source Model Image"),
                    colormap=colormap, vmax=source_vmax,
                    lines=image_plane_lines, line_colors=image_plane_line_colors)
     else:
         axes_flat[1].axis("off")
 
-
     # Panel 2: Source Plane (No Zoom) (same as subplot_fit panel 12)
     _plot_source_plane(fit, axes_flat[2], final_plane_index, zoom_to_brightest=False,
-                       colormap=colormap, title="Source Plane (No Zoom)",
+                       colormap=colormap, title=_pf("Source Plane (No Zoom)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors,
                        vmax=source_vmax)
 
     # Panel 3: Lens Image (log10)
-    plot_array(array=lens_image, ax=axes_flat[3], title="Lens Image",
+    plot_array(array=lens_image, ax=axes_flat[3], title=_pf("Lens Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                colormap=colormap, use_log10=True)
 
     # Panel 4: Convergence (log10)
     try:
         convergence = tracer.convergence_2d_from(grid=grid)
-        plot_array(array=convergence, ax=axes_flat[4], title="Convergence",
+        plot_array(array=convergence, ax=axes_flat[4], title=_pf("Convergence"),
                    colormap=colormap, use_log10=True)
     except Exception:
         axes_flat[4].axis("off")
@@ -744,21 +760,21 @@ def subplot_tracer_from_fit(
     # Panel 5: Potential (log10)
     try:
         potential = tracer.potential_2d_from(grid=grid)
-        plot_array(array=potential, ax=axes_flat[5], title="Potential",
+        plot_array(array=potential, ax=axes_flat[5], title=_pf("Potential"),
                    colormap=colormap, use_log10=True)
     except Exception:
         axes_flat[5].axis("off")
 
     # Panel 6: Deflections Y
-    plot_array(array=deflections_y, ax=axes_flat[6], title="Deflections Y",
+    plot_array(array=deflections_y, ax=axes_flat[6], title=_pf("Deflections Y"),
                colormap=colormap)
 
     # Panel 7: Deflections X
-    plot_array(array=deflections_x, ax=axes_flat[7], title="Deflections X",
+    plot_array(array=deflections_x, ax=axes_flat[7], title=_pf("Deflections X"),
                colormap=colormap)
 
     # Panel 8: Magnification
-    plot_array(array=magnification, ax=axes_flat[8], title="Magnification",
+    plot_array(array=magnification, ax=axes_flat[8], title=_pf("Magnification"),
                colormap=colormap)
 
     tight_layout()
@@ -770,6 +786,7 @@ def subplot_fit_combined(
     output_path: Optional[str] = None,
     output_format: str = None,
     colormap: Optional[str] = None,
+    title_prefix: str = None,
 ):
     """
     Produce a combined multi-row subplot for a list of `FitImaging` objects.
@@ -809,39 +826,41 @@ def subplot_fit_combined(
 
     final_plane_index = len(fit_list[0].tracer.planes) - 1
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     for row, fit in enumerate(fit_list):
         row_axes = all_axes[row]
 
-        plot_array(array=fit.data, ax=row_axes[0], title="Data", colormap=colormap)
+        plot_array(array=fit.data, ax=row_axes[0], title=_pf("Data"), colormap=colormap)
 
         try:
             subtracted = fit.subtracted_images_of_planes_list[1]
-            plot_array(array=subtracted, ax=row_axes[1], title="Subtracted Image",
+            plot_array(array=subtracted, ax=row_axes[1], title=_pf("Subtracted Image"),
                        colormap=colormap)
         except (IndexError, AttributeError):
             row_axes[1].axis("off")
 
         try:
             lens_model = fit.model_images_of_planes_list[0]
-            plot_array(array=lens_model, ax=row_axes[2], title="Lens Model Image",
+            plot_array(array=lens_model, ax=row_axes[2], title=_pf("Lens Model Image"),
                        colormap=colormap)
         except (IndexError, AttributeError):
             row_axes[2].axis("off")
 
         try:
             source_model = fit.model_images_of_planes_list[final_plane_index]
-            plot_array(array=source_model, ax=row_axes[3], title="Source Model Image",
+            plot_array(array=source_model, ax=row_axes[3], title=_pf("Source Model Image"),
                        colormap=colormap)
         except (IndexError, AttributeError):
             row_axes[3].axis("off")
 
         try:
-            _plot_source_plane(fit, row_axes[4], final_plane_index, colormap=colormap)
+            _plot_source_plane(fit, row_axes[4], final_plane_index, colormap=colormap,
+                               title=_pf(f"Source Plane {final_plane_index}"))
         except Exception:
             row_axes[4].axis("off")
 
         plot_array(array=fit.normalized_residual_map, ax=row_axes[5],
-                   title="Normalized Residual Map", colormap=colormap, cb_unit=r"$\sigma$")
+                   title=_pf("Normalized Residual Map"), colormap=colormap, cb_unit=r"$\sigma$")
 
     tight_layout()
     save_figure(fig, path=output_path, filename="fit_combined", format=output_format)
@@ -852,6 +871,7 @@ def subplot_fit_combined_log10(
     output_path: Optional[str] = None,
     output_format: str = None,
     colormap: Optional[str] = None,
+    title_prefix: str = None,
 ):
     """
     Produce a combined log10 multi-row subplot for a list of `FitImaging` objects.
@@ -883,41 +903,43 @@ def subplot_fit_combined_log10(
 
     final_plane_index = len(fit_list[0].tracer.planes) - 1
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     for row, fit in enumerate(fit_list):
         row_axes = all_axes[row]
 
-        plot_array(array=fit.data, ax=row_axes[0], title="Data", colormap=colormap,
+        plot_array(array=fit.data, ax=row_axes[0], title=_pf("Data"), colormap=colormap,
                    use_log10=True)
 
         try:
             subtracted = fit.subtracted_images_of_planes_list[1]
-            plot_array(array=subtracted, ax=row_axes[1], title="Subtracted Image",
+            plot_array(array=subtracted, ax=row_axes[1], title=_pf("Subtracted Image"),
                        colormap=colormap, use_log10=True)
         except (IndexError, AttributeError):
             row_axes[1].axis("off")
 
         try:
             lens_model = fit.model_images_of_planes_list[0]
-            plot_array(array=lens_model, ax=row_axes[2], title="Lens Model Image",
+            plot_array(array=lens_model, ax=row_axes[2], title=_pf("Lens Model Image"),
                        colormap=colormap, use_log10=True)
         except (IndexError, AttributeError):
             row_axes[2].axis("off")
 
         try:
             source_model = fit.model_images_of_planes_list[final_plane_index]
-            plot_array(array=source_model, ax=row_axes[3], title="Source Model Image",
+            plot_array(array=source_model, ax=row_axes[3], title=_pf("Source Model Image"),
                        colormap=colormap, use_log10=True)
         except (IndexError, AttributeError):
             row_axes[3].axis("off")
 
         try:
             _plot_source_plane(fit, row_axes[4], final_plane_index, colormap=colormap,
-                               use_log10=True)
+                               use_log10=True,
+                               title=_pf(f"Source Plane {final_plane_index}"))
         except Exception:
             row_axes[4].axis("off")
 
         plot_array(array=fit.normalized_residual_map, ax=row_axes[5],
-                   title="Normalized Residual Map", colormap=colormap, cb_unit=r"$\sigma$")
+                   title=_pf("Normalized Residual Map"), colormap=colormap, cb_unit=r"$\sigma$")
 
     tight_layout()
     save_figure(fig, path=output_path, filename="fit_combined_log10", format=output_format)
diff --git a/autolens/interferometer/model/plotter.py b/autolens/interferometer/model/plotter.py
index 4d7b7c4ee..53caf0572 100644
--- a/autolens/interferometer/model/plotter.py
+++ b/autolens/interferometer/model/plotter.py
@@ -73,6 +73,7 @@ def should_plot(name):
                 fit, output_path=output_path, output_format=fmt,
                 image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                 source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                title_prefix=self.title_prefix,
             )
 
         if plot_setting(section="tracer", name="subplot_tracer"):
@@ -80,12 +81,14 @@ def should_plot(name):
                 fit, output_path=output_path, output_format=fmt,
                 image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
                 source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                title_prefix=self.title_prefix,
             )
 
         if should_plot("subplot_fit_dirty_images") or quick_update:
             subplot_fit_dirty_images(
                 fit, output_path=output_path, output_format=fmt,
                 image_plane_lines=ip_lines, image_plane_line_colors=ip_colors,
+                title_prefix=self.title_prefix,
             )
 
         if quick_update:
@@ -95,6 +98,7 @@ def should_plot(name):
             subplot_fit_real_space(
                 fit, output_path=output_path, output_format=fmt,
                 source_plane_lines=sp_lines, source_plane_line_colors=sp_colors,
+                title_prefix=self.title_prefix,
             )
 
         if should_plot("fits_galaxy_images"):
diff --git a/autolens/interferometer/plot/fit_interferometer_plots.py b/autolens/interferometer/plot/fit_interferometer_plots.py
index 7a673ba23..8b33ffe9d 100644
--- a/autolens/interferometer/plot/fit_interferometer_plots.py
+++ b/autolens/interferometer/plot/fit_interferometer_plots.py
@@ -142,6 +142,7 @@ def subplot_fit(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 12-panel subplot summarising an interferometer fit.
@@ -183,6 +184,7 @@ def subplot_fit(
             _compute_critical_curve_lines(tracer, _cc_grid)
         )
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 4, figsize=conf_subplot_figsize(3, 4))
     axes_flat = list(axes.flatten())
 
@@ -191,32 +193,32 @@ def subplot_fit(
         y=np.real(fit.residual_map),
         x=fit.dataset.uv_distances / 10 ** 3.0,
         ax=axes_flat[0],
-        title="Amplitudes vs UV-Distance",
+        title=_pf("Amplitudes vs UV-Distance"),
         xtick_suffix='"',
         ytick_suffix="Jy",
         plot_axis_type="scatter",
     )
 
-    plot_array(array=fit.dirty_image, ax=axes_flat[1], title="Dirty Image",
+    plot_array(array=fit.dirty_image, ax=axes_flat[1], title=_pf("Dirty Image"),
                colormap=colormap)
     plot_array(array=fit.dirty_signal_to_noise_map, ax=axes_flat[2],
-               title="Dirty Signal-To-Noise Map", colormap=colormap)
+               title=_pf("Dirty Signal-To-Noise Map"), colormap=colormap)
 
     # Panel 3 (4th): dirty model image with critical curves
-    plot_array(array=fit.dirty_model_image, ax=axes_flat[3], title="Dirty Model Image",
+    plot_array(array=fit.dirty_model_image, ax=axes_flat[3], title=_pf("Dirty Model Image"),
                colormap=colormap, lines=image_plane_lines,
                line_colors=image_plane_line_colors)
 
     # Panel 4: dirty residual map
     plot_array(array=fit.dirty_residual_map, ax=axes_flat[4],
-               title="Dirty Residual Map", colormap=colormap)
+               title=_pf("Dirty Residual Map"), colormap=colormap)
 
     # Panel 5: normalized residual vs UV-distances (real)
     plot_yx(
         y=np.real(fit.normalized_residual_map),
         x=fit.dataset.uv_distances / 10 ** 3.0,
         ax=axes_flat[5],
-        title="Normalized Residual Map (Real)",
+        title=_pf("Normalized Residual Map (Real)"),
         xtick_suffix='"',
         ytick_suffix=r"$\sigma$",
         plot_axis_type="scatter",
@@ -227,7 +229,7 @@ def subplot_fit(
         y=np.imag(fit.normalized_residual_map),
         x=fit.dataset.uv_distances / 10 ** 3.0,
         ax=axes_flat[6],
-        title="Normalized Residual Map (Imag)",
+        title=_pf("Normalized Residual Map (Imag)"),
         xtick_suffix='"',
         ytick_suffix=r"$\sigma$",
         plot_axis_type="scatter",
@@ -236,30 +238,30 @@ def subplot_fit(
     # Panel 7 (8th): source plane zoomed with caustics
     _plot_source_plane(fit, axes_flat[7], final_plane_index,
                        zoom_to_brightest=True, colormap=colormap,
-                       title="Source Plane (Zoomed)",
+                       title=_pf("Source Plane (Zoomed)"),
                        lines=source_plane_lines,
                        line_colors=source_plane_line_colors)
 
     plot_array(array=fit.dirty_normalized_residual_map, ax=axes_flat[8],
-               title="Dirty Normalized Residual Map", colormap=colormap, cb_unit=r"$\sigma$")
+               title=_pf("Dirty Normalized Residual Map"), colormap=colormap, cb_unit=r"$\sigma$")
 
     # Panel 9: clipped to ±1σ
     plot_array(
         fit.dirty_normalized_residual_map,
         ax=axes_flat[9],
-        title=r"Normalized Residual Map $1\sigma$",
+        title=_pf(r"Normalized Residual Map $1\sigma$"),
         colormap=colormap,
         vmin=-1.0, vmax=1.0,
         cb_unit=r"$\sigma$",
     )
 
     plot_array(array=fit.dirty_chi_squared_map, ax=axes_flat[10],
-               title="Dirty Chi-Squared Map", colormap=colormap, cb_unit=r"$\chi^2$")
+               title=_pf("Dirty Chi-Squared Map"), colormap=colormap, cb_unit=r"$\chi^2$")
 
     # Panel 11 (12th): source plane not zoomed with caustics
     _plot_source_plane(fit, axes_flat[11], final_plane_index,
                        zoom_to_brightest=False, colormap=colormap,
-                       title="Source Plane (No Zoom)",
+                       title=_pf("Source Plane (No Zoom)"),
                        lines=source_plane_lines,
                        line_colors=source_plane_line_colors)
 
@@ -275,6 +277,7 @@ def subplot_fit_dirty_images(
     use_log10: bool = False,
     image_plane_lines=None,
     image_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 2×3 subplot of dirty-image diagnostics for an interferometer fit.
@@ -303,22 +306,23 @@ def subplot_fit_dirty_images(
             _compute_critical_curve_lines(tracer, _cc_grid)
         )
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(2, 3, figsize=conf_subplot_figsize(2, 3))
     axes_flat = list(axes.flatten())
 
-    plot_array(array=fit.dirty_image, ax=axes_flat[0], title="Dirty Image",
+    plot_array(array=fit.dirty_image, ax=axes_flat[0], title=_pf("Dirty Image"),
                colormap=colormap, use_log10=use_log10)
     plot_array(array=fit.dirty_signal_to_noise_map, ax=axes_flat[1],
-               title="Dirty Signal-To-Noise Map", colormap=colormap)
+               title=_pf("Dirty Signal-To-Noise Map"), colormap=colormap)
     plot_array(array=fit.dirty_model_image, ax=axes_flat[2],
-               title="Dirty Model Image", colormap=colormap, use_log10=use_log10,
+               title=_pf("Dirty Model Image"), colormap=colormap, use_log10=use_log10,
                lines=image_plane_lines, line_colors=image_plane_line_colors)
     plot_array(array=fit.dirty_residual_map, ax=axes_flat[3],
-               title="Dirty Residual Map", colormap=colormap)
+               title=_pf("Dirty Residual Map"), colormap=colormap)
     plot_array(array=fit.dirty_normalized_residual_map, ax=axes_flat[4],
-               title="Dirty Normalized Residual Map", colormap=colormap, cb_unit=r"$\sigma$")
+               title=_pf("Dirty Normalized Residual Map"), colormap=colormap, cb_unit=r"$\sigma$")
     plot_array(array=fit.dirty_chi_squared_map, ax=axes_flat[5],
-               title="Dirty Chi-Squared Map", colormap=colormap, cb_unit=r"$\chi^2$")
+               title=_pf("Dirty Chi-Squared Map"), colormap=colormap, cb_unit=r"$\chi^2$")
 
     tight_layout()
     save_figure(fig, path=output_path, filename="fit_dirty_images", format=output_format)
@@ -331,6 +335,7 @@ def subplot_fit_real_space(
     colormap: Optional[str] = None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a real-space subplot for an interferometer fit.
@@ -363,23 +368,24 @@ def subplot_fit_real_space(
     fig, axes = subplots(1, 2, figsize=conf_subplot_figsize(1, 2))
     axes_flat = list(axes.flatten())
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     if fit.inversion is None:
         # Parametric source: image-plane model image + source-plane image
         grid = fit.dataset.real_space_mask.derive_grid.all_false
         image = tracer.image_2d_from(grid=grid)
-        plot_array(array=image, ax=axes_flat[0], title="Image", colormap=colormap)
+        plot_array(array=image, ax=axes_flat[0], title=_pf("Image"), colormap=colormap)
 
         _plot_source_plane(fit, axes_flat[1], final_plane_index,
                            zoom_to_brightest=True, colormap=colormap,
-                           title="Source Plane (Zoomed)",
+                           title=_pf("Source Plane (Zoomed)"),
                            lines=source_plane_lines, line_colors=source_plane_line_colors)
     else:
         # Pixelized source: dirty model image + source reconstruction
         plot_array(array=fit.dirty_model_image, ax=axes_flat[0],
-                   title="Reconstructed Image", colormap=colormap)
+                   title=_pf("Reconstructed Image"), colormap=colormap)
         _plot_source_plane(fit, axes_flat[1], final_plane_index,
                            zoom_to_brightest=True, colormap=colormap,
-                           title="Source Reconstruction",
+                           title=_pf("Source Reconstruction"),
                            lines=source_plane_lines, line_colors=source_plane_line_colors)
 
     tight_layout()
@@ -395,6 +401,7 @@ def subplot_tracer_from_fit(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a 9-panel tracer subplot derived from a `FitInterferometer` object.
@@ -445,11 +452,12 @@ def subplot_tracer_from_fit(
 
     magnification = LensCalc.from_mass_obj(tracer).magnification_2d_from(grid=grid)
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 3, figsize=conf_subplot_figsize(3, 3))
     axes_flat = list(axes.flatten())
 
     # Panel 0: Dirty Model Image
-    plot_array(array=fit.dirty_model_image, ax=axes_flat[0], title="Dirty Model Image",
+    plot_array(array=fit.dirty_model_image, ax=axes_flat[0], title=_pf("Dirty Model Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                colormap=colormap)
 
@@ -468,7 +476,7 @@ def subplot_tracer_from_fit(
     except Exception:
         source_model_img = None
     if source_model_img is not None:
-        plot_array(array=source_model_img, ax=axes_flat[1], title="Source Model Image",
+        plot_array(array=source_model_img, ax=axes_flat[1], title=_pf("Source Model Image"),
                    colormap=colormap,
                    lines=image_plane_lines, line_colors=image_plane_line_colors)
     else:
@@ -476,18 +484,18 @@ def subplot_tracer_from_fit(
 
     # Panel 2: Source Plane (No Zoom) with caustics
     _plot_source_plane(fit, axes_flat[2], final_plane_index, zoom_to_brightest=False,
-                       colormap=colormap, title="Source Plane (No Zoom)",
+                       colormap=colormap, title=_pf("Source Plane (No Zoom)"),
                        lines=source_plane_lines, line_colors=source_plane_line_colors)
 
     # Panel 3: Lens Image (log10)
-    plot_array(array=lens_image, ax=axes_flat[3], title="Lens Image",
+    plot_array(array=lens_image, ax=axes_flat[3], title=_pf("Lens Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                colormap=colormap, use_log10=True)
 
     # Panel 4: Convergence (log10)
     try:
         convergence = tracer.convergence_2d_from(grid=grid)
-        plot_array(array=convergence, ax=axes_flat[4], title="Convergence",
+        plot_array(array=convergence, ax=axes_flat[4], title=_pf("Convergence"),
                    colormap=colormap, use_log10=True)
     except Exception:
         axes_flat[4].axis("off")
@@ -495,21 +503,21 @@ def subplot_tracer_from_fit(
     # Panel 5: Potential (log10)
     try:
         potential = tracer.potential_2d_from(grid=grid)
-        plot_array(array=potential, ax=axes_flat[5], title="Potential",
+        plot_array(array=potential, ax=axes_flat[5], title=_pf("Potential"),
                    colormap=colormap, use_log10=True)
     except Exception:
         axes_flat[5].axis("off")
 
     # Panel 6: Deflections Y
-    plot_array(array=deflections_y, ax=axes_flat[6], title="Deflections Y",
+    plot_array(array=deflections_y, ax=axes_flat[6], title=_pf("Deflections Y"),
                colormap=colormap)
 
     # Panel 7: Deflections X
-    plot_array(array=deflections_x, ax=axes_flat[7], title="Deflections X",
+    plot_array(array=deflections_x, ax=axes_flat[7], title=_pf("Deflections X"),
                colormap=colormap)
 
     # Panel 8: Magnification
-    plot_array(array=magnification, ax=axes_flat[8], title="Magnification",
+    plot_array(array=magnification, ax=axes_flat[8], title=_pf("Magnification"),
                colormap=colormap)
 
     tight_layout()
diff --git a/autolens/lens/plot/tracer_plots.py b/autolens/lens/plot/tracer_plots.py
index 2cf6bd78d..79df60add 100644
--- a/autolens/lens/plot/tracer_plots.py
+++ b/autolens/lens/plot/tracer_plots.py
@@ -112,6 +112,7 @@ def subplot_tracer(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """Multi-panel subplot of the tracer: image, source images, and mass quantities.
 
@@ -160,29 +161,30 @@ def subplot_tracer(
 
     magnification = LensCalc.from_mass_obj(tracer).magnification_2d_from(grid=grid)
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     fig, axes = subplots(3, 3, figsize=conf_subplot_figsize(3, 3))
     axes_flat = list(axes.flatten())
 
-    plot_array(array=image, ax=axes_flat[0], title="Model Image",
+    plot_array(array=image, ax=axes_flat[0], title=_pf("Model Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                positions=pos_list, colormap=colormap, use_log10=use_log10)
-    plot_array(array=source_image, ax=axes_flat[1], title="Source Model Image",
+    plot_array(array=source_image, ax=axes_flat[1], title=_pf("Source Model Image"),
                colormap=colormap, use_log10=use_log10, vmax=source_vmax)
-    plot_array(array=source_image, ax=axes_flat[2], title="Source Plane (No Zoom)",
+    plot_array(array=source_image, ax=axes_flat[2], title=_pf("Source Plane (No Zoom)"),
                lines=source_plane_lines, line_colors=source_plane_line_colors,
                colormap=colormap, use_log10=use_log10)
-    plot_array(array=lens_image, ax=axes_flat[3], title="Lens Image",
+    plot_array(array=lens_image, ax=axes_flat[3], title=_pf("Lens Image"),
                lines=image_plane_lines, line_colors=image_plane_line_colors,
                colormap=colormap, use_log10=True)
-    plot_array(array=convergence, ax=axes_flat[4], title="Convergence",
+    plot_array(array=convergence, ax=axes_flat[4], title=_pf("Convergence"),
                colormap=colormap, use_log10=True)
-    plot_array(array=potential, ax=axes_flat[5], title="Potential",
+    plot_array(array=potential, ax=axes_flat[5], title=_pf("Potential"),
                colormap=colormap, use_log10=True)
-    plot_array(array=deflections_y, ax=axes_flat[6], title="Deflections Y",
+    plot_array(array=deflections_y, ax=axes_flat[6], title=_pf("Deflections Y"),
                lines=image_plane_lines, line_colors=image_plane_line_colors, colormap=colormap)
-    plot_array(array=deflections_x, ax=axes_flat[7], title="Deflections X",
+    plot_array(array=deflections_x, ax=axes_flat[7], title=_pf("Deflections X"),
                lines=image_plane_lines, line_colors=image_plane_line_colors, colormap=colormap)
-    plot_array(array=magnification, ax=axes_flat[8], title="Magnification",
+    plot_array(array=magnification, ax=axes_flat[8], title=_pf("Magnification"),
                lines=image_plane_lines, line_colors=image_plane_line_colors, colormap=colormap)
 
     hide_unused_axes(axes_flat)
@@ -197,6 +199,7 @@ def subplot_lensed_images(
     output_format: str = None,
     colormap: Optional[str] = None,
     use_log10: bool = False,
+    title_prefix: str = None,
 ):
     """
     Produce a subplot with one panel per tracer plane showing each plane's image.
@@ -229,13 +232,14 @@ def subplot_lensed_images(
     fig, axes = subplots(1, n, figsize=conf_subplot_figsize(1, n))
     axes_flat = [axes] if n == 1 else list(np.array(axes).flatten())
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     for plane_index in range(n):
         galaxies = ag.Galaxies(galaxies=tracer.planes[plane_index])
         image = galaxies.image_2d_from(grid=traced_grids[plane_index])
         plot_array(
             array=image,
             ax=axes_flat[plane_index],
-            title=f"Image Of Plane {plane_index}",
+            title=_pf(f"Image Of Plane {plane_index}"),
             colormap=colormap,
             use_log10=use_log10,
         )
@@ -251,6 +255,7 @@ def subplot_galaxies_images(
     output_format: str = None,
     colormap: Optional[str] = None,
     use_log10: bool = False,
+    title_prefix: str = None,
 ):
     """
     Produce a subplot showing per-galaxy images for every plane in the tracer.
@@ -291,12 +296,13 @@ def subplot_galaxies_images(
 
     idx = 0
 
+    _pf = (lambda t: f"{title_prefix.rstrip()} {t}") if title_prefix else (lambda t: t)
     lens_galaxies = ag.Galaxies(galaxies=tracer.planes[0])
     lens_image = lens_galaxies.image_2d_from(grid=traced_grids[0])
     plot_array(
         array=lens_image,
         ax=axes_flat[idx],
-        title="Image Of Plane 0",
+        title=_pf("Image Of Plane 0"),
         colormap=colormap,
         use_log10=use_log10,
     )
@@ -311,7 +317,7 @@ def subplot_galaxies_images(
             plot_array(
                 array=image,
                 ax=axes_flat[idx],
-                title=f"Image Of Plane {plane_index}",
+                title=_pf(f"Image Of Plane {plane_index}"),
                 colormap=colormap,
                 use_log10=use_log10,
             )
@@ -321,7 +327,7 @@ def subplot_galaxies_images(
             plot_array(
                 array=image,
                 ax=axes_flat[idx],
-                title=f"Plane Image Of Plane {plane_index}",
+                title=_pf(f"Plane Image Of Plane {plane_index}"),
                 colormap=colormap,
                 use_log10=use_log10,
             )
diff --git a/autolens/point/model/plotter.py b/autolens/point/model/plotter.py
index add77a389..74fe8d941 100644
--- a/autolens/point/model/plotter.py
+++ b/autolens/point/model/plotter.py
@@ -29,7 +29,8 @@ def should_plot(name):
         fmt = self.fmt
 
         if should_plot("subplot_dataset"):
-            subplot_dataset(dataset, output_path=output_path, output_format=fmt)
+            subplot_dataset(dataset, output_path=output_path, output_format=fmt,
+                            title_prefix=self.title_prefix)
 
     def fit_point(
         self,
@@ -84,6 +85,7 @@ def should_plot(name):
                 image_plane_line_colors=ip_colors,
                 source_plane_lines=sp_lines,
                 source_plane_line_colors=sp_colors,
+                title_prefix=self.title_prefix,
             )
 
         if quick_update:
diff --git a/autolens/point/plot/fit_point_plots.py b/autolens/point/plot/fit_point_plots.py
index 4201b8175..b9f190837 100644
--- a/autolens/point/plot/fit_point_plots.py
+++ b/autolens/point/plot/fit_point_plots.py
@@ -12,6 +12,7 @@ def subplot_fit(
     image_plane_line_colors=None,
     source_plane_lines=None,
     source_plane_line_colors=None,
+    title_prefix: str = None,
 ):
     """
     Produce a subplot summarising a `FitPointDataset`.
@@ -56,10 +57,11 @@ def subplot_fit(
         else fit.positions.model_data
     )
 
+    _prefix = f"{title_prefix.rstrip()} " if title_prefix else ""
     plot_grid(
         grid=obs_grid,
         ax=axes_flat[0],
-        title=f"{fit.dataset.name} Fit Positions",
+        title=f"{_prefix}{fit.dataset.name} Fit Positions",
         output_path=None,
         output_filename=None,
         output_format=output_format,
@@ -74,7 +76,7 @@ def subplot_fit(
             y=y,
             x=x,
             ax=axes_flat[1],
-            title=f"{fit.dataset.name} Fit Fluxes",
+            title=f"{_prefix}{fit.dataset.name} Fit Fluxes",
             output_path=None,
             output_filename="fit_point_fluxes",
             output_format=output_format,
diff --git a/autolens/point/plot/point_dataset_plots.py b/autolens/point/plot/point_dataset_plots.py
index 4fe6c1794..7ca54aa0f 100644
--- a/autolens/point/plot/point_dataset_plots.py
+++ b/autolens/point/plot/point_dataset_plots.py
@@ -8,6 +8,7 @@ def subplot_dataset(
     dataset,
     output_path: Optional[str] = None,
     output_format: str = None,
+    title_prefix: str = None,
 ):
     """
     Produce a subplot visualising a `PointDataset`.
@@ -46,10 +47,11 @@ def subplot_dataset(
         else dataset.positions
     )
 
+    _prefix = f"{title_prefix.rstrip()} " if title_prefix else ""
     plot_grid(
         grid=grid,
         ax=axes_flat[0],
-        title=f"{dataset.name} Positions",
+        title=f"{_prefix}{dataset.name} Positions",
         output_path=None,
         output_filename=None,
         output_format=output_format,
@@ -62,7 +64,7 @@ def subplot_dataset(
             y=y,
             x=x,
             ax=axes_flat[1],
-            title=f"{dataset.name} Fluxes",
+            title=f"{_prefix}{dataset.name} Fluxes",
             output_path=None,
             output_filename="point_dataset_fluxes",
             output_format=output_format,